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Dynamic Response of Heat Exchangers 
Having Internal Heat Sources—Part Ill 


VEDAT S. ARPACI This paper is an extension of two of the authors’ previous papers |1, 2|* and considers 

the dynamic response of the temperature of the heat-transfer surface and the fluid-surface 

Instructor, Mechanical Engineering temperature difference in a heat exchanger having internal heat sources. A general 
Faculty, Istanbul pag ohn yn solution is presented for a transient resulting from an arbitrary time rate of change of 
ae ne heat generation starting from an arbitrary initial condition. The results are obtained 


by a direct mathematical attack on the governing differential equations. Asa special 
JOHN A. CLARK case, the transient resulting from a linear time rate of change of inter) al heat generation 
Professor of Mechanical Engineering, is given, and this ts compared with both experimental and theoretical results correspond 
University of Michigan, Ann Arbor, ing to a step change in heat generation. Satisfactory comparison is obtained. The re 
Mich. Assoc. Mem. ASME sults for this special case are presented such that the surface temperature and fluid 
surface temperature difference may be computed as a function of space and time using 
four dimensionless, mathematical functions which are given in graphical form. It is 
shown haw the case of the linear time rate of change of internal heat generation may be 
used to compute the dynamic response of this type of heat exchanger for any transient 
resulting from excursions in heat-generation rate which may be approximated piece 
wise by a series of linear functions of time. Heat exchangers to which these solutions 
apply include the heterogeneous nuclear reactor. 


Introduction mathematical solution to the response of the temperature of the 

a reactor coolant. As a consequence of this study, the solution for 
ROBLEMs ASSOCIATED with the dynamic behavior of the response of this type of heat exchanger to several different 
systems have application in the field of automatic control which is kinds of induced transients was obtained and has resulted in this 
becoming of increasing technical importance. It was in this con-  geries of papers. These results apply, of course, to other types of 
nection that, during the summer of 1954, the authors became in- heat exchangers having internal heat sources, including an elec- 
terested in the dynamic behavior of heat exchangers, particularly trical heater and a chemical reactor in which a chemical reaction 
those having internal heat sources. At that time it was necessary 


to determine the time response of a fluid circuit contaiming a nu- 


occurs within the solid walls. Apart from the use in the control 


of the systems, knowledge of transient behavior is important in 
‘ar reactor, piping, valves, and heat exchangers in which the design where thermal stresses are to be considered or when tem- 


disturbance was caused by power transients in the reactor The perature-time information is required from a device such as a heat 


reactor was characterized as a heat exchanger having internal — exehanger 
heat sources One approach included an atte mpt to obtain a This is the third in a series of papers which has been planned 
on the general subject of the dynamic response of heat exchangers 


Present address: Department of Mecha Engineering, M.I.T., having time-variant internal heat sources. The first two papers 
Cambridee. Mass gp ; a 
- & , ; s 1, 2] were presented in June and August, 1957. 
? Numbers in brackets designate Referé it end of paper ; 
Contributed by the Heat Transfer Divi nd presented at the 


Semi-Annual Meeting, Detroit, Mich., Ju 5-19, 1958. of Tue program such problems as these on digital or analog computing 


It is recognized that often it is convenient and practical to 


(AMERICAN SOCIETY OF MECHANICAL ENGINEER machines, particularly if a great many different cases are to be 
Statements and opinions advanced in papers are to be j ; 

atemen and opinion papers are » be studied However, the programming is by no means trivial and 

as individual expressions of their ithors and not those 

Manuscript received at ASME Headquarters, April 7, 

No. 58—SA-39 solutions. Furthermore, for those having access to machines, 


machine time is not always available or possible to those seeking 





Nomenclature 


= area of the fluid in contact with wall, sq ft fi, fe = funetions defined in Appendix 2 
flow area of fluid, sq ft 1a) -- -y Fig = see Appendix 2 
cross-sectional area of wal | 13> ,Gig = see Appendix 2 
defined by equation (3a eg coefficient of heat transfer between the wall and 
defined by equation (3b), hr the fluid, Btu/hr-ft*-deg I 
defined by equation (3a : enthalpy, Btu/lby 
defined by equation (3 it = Bessel function of first kind, zeroth order 
a, +a; = (M + 1)/K, hi k = the kth term of the series 
Qo, = M/K2, (hr ‘ (h/pe,)A/V, (hour)~', for fluid 
circumference (wetted perin t w = A/(pe,)y (A/V), (hour , for wall 
specific heat of the fluid, Btu /IlI mass, lby 
specific heat of the wall, Btu/! g 203, heat capacity ratio, pc, \ w/ Pe 3 
outside diameter of the tub I = pressure, psia 
inside diameter of the tube, { Laplace variable, (hour 
internal energy, Btu/Iby, « } . = rate of work done by sys 


base of natural logarithms, dimensionless (Continued on next page) 
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these sol serve as helpful check solutions, and for those 


itions may 


who do not, the solutions presented here may be a useful tool 


iny event, 


with a small expenditure of time and at little ex- 


pense 


i human computer using a desk calculator may run out a 


number ises in order to demonstrate the effect of a chosen 


on the response, thus giving a ‘feel’ for this type of 


variable 


dynamic system. Actually, machine computation was employed 


in these solutions as most of the functions which are given were 


obtained on an IBM 650 digital computer. Only the function 


¥ie(s, g*) was reduced by hand computation 


Literature Survey 


Since the first two papers a several new 


To 


the literat ire previously reviewed 


writing of the 


papers on this top have come to the authors’ attention 


ivold unnecessary repetition 


will not be included in the body of this paper, but for convenience 


the References are listed at the end Comment will be made on 


the new material 
The 
! 


heat exchangers is partially 


sis of the on of counterflow and parallel flow 
11) and solved by 
Rizika [11] 
has computed and presented as graphs two infinite series, one of 


y.** 


\ numerical method for calculating heat- 


— 
ria y 
treated by Rizika 


(arslaw and Jaeger [10) for zero wall-heat capacity 


which is employe d in this paper as part of the solution to 
the present proble n 
exchanger dynamics is given by Dusinberre [12] who presents 
Detailed 


\ problem allied with that of the present 


explicit iteration formulas and computation guides 


examples are given 
paper is published by Fritz [13] in which the response of a nuclear- 


reactor coolant owing to time-varying flow rates is given for a 


lumped-temperature system. Thermal response of a homogene- 
reactor has been treated by Gallagher [14] 
ind Dietrich 


ind coolant te mperatures in @ heterogeneous nuclear reactor 


ous type nuclear 


Silvers 15) have computed the response ol the sur- 
lace 
resulting from a step increase in power using the method of finite 
The if 


sponding to the exact mathematical solution given in Figs. 7 and 8 


differences results bear exact similarity to those corre- 


ol reterence |2 


An extension of the problem studied by Rizika [9] and Carslaw 
and Jaeger [10] has been made by Dabora [16] in which the 
effect of heat loss to the surroundings from a regenerative-type 
The 
obtained by ignoring the time rate of change of fluid tempera- 
ture, similar to that of Rizika [11] 
flow heat exchanger (rotary regenerator) has been computed by 


Nine 


tables of numerical results corresponding to a fairly wide range 


(nonperiodic) heat exchanger is considered solution is 


The effectiveness of a periodic- 


Lambertson [17] using the finite-difference technique. 


of operating conditions are given. The effect of heat transfer be- 
tween a gas and the walls and other solid material in a container 
tank on the charging and blowdown processes is examined by 
Reynolds and Kays [18]. 


cases of constant-mass flow and blowdown through a critical-flow 


The authors have considered the two 


nozzle. Analytical solutions assuming a lumped-parameter sys- 
tem are presented which agree favorably with experimental re- 


Cima and London [19] have studied experimentally the 


sults [ 
transient response of a particular two-fluid counterflow heat ex- 
( hanger using an electrical analog computer having four discrete 
lumped circuits corresponding to equal increments of flow length 
a gas-turbine re- 


in the exchanger. The heat exchanger studied is 


generator of equal hot and cold-side fluid-<« pacity rates and 


equal Dimensionless results are 


thermal resistances given for 
the response of the outlet cold-gas temperature to step changes in 
the inlet hot-gas temperature and for step changes in the flow 
rate of the cold 


estimation of the response of the cold-gas outlet temperature to 


and hot—gas. A procedure is discussed for the 
an arbitrary disturbance in the hot-gas inlet temperature and re- 
sults are given for a particular case. 

Other work concerning the transient 
fluids is that of Rea and Ablow [20]. Their sy 
an insulated duct through which heated air was passed for condi- 


recent flow of heated 


stem consisted of 


tions approximating a sudden and oscillating change in both flow 
and temperature of the inlet air. Measurements are reported of 
the response of the wall temperature and the inlet and outlet air 
temperatures and these are compared with a finite-difference 


calculation. Siegel [21] has reported on the transient process in 


free convection from a semi-infinite vertical flat plate, initially at 





Nomenclature 


yy initial power input, Btu/ft® hr 
Ap,’ change in power input, Btu/ft*® hr 
q/A heat flux, Btu/hr-ft? 
q/A initial heat flux, Btu/hr-ft? 
heat flow rate, Btu/hr 
outside radius of the tube, ft 
inside radius of the tube, ft 
K/M) x/u, 


K(2/u), dimensionless 


dimensionless 


temperature of the fluid, deg F 
(x, 7) 

state component of fluid temperature, deg 
transient component of fluid temperature, deg F 
temperature of the fluid entering coolant channel, 

deg | 

velocity of the fluid, fps 
Ibu 
volume of the fluid, ft? 
volume of the wall, ft® 
fluid flow rate 


volume, ft 


aper ifie 


Iby /hr 

axial distance, ft 

, lbw /ft 

transient in 
Btu /hr-ft® 

wall temperature, deg F, @(x, r 
6,(x, 0) + Ox, r) 


density 


volumetric heat-generation rate, 
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6,(x, 0) = steady-state component of wall temperature, deg 


transient component of wall temperature, deg F 
time, hr 
rT — z/u, hr 
characteristic time for linear transient, hr, see 
equations (11) 
see equation (10¢ 
see equation (10f) 
6 t, deg F 
Appendix 2, dimensionless 
Appendix 2, dimensionless 
) 


Appendix 2, dimensionless | 


2, dimensionless 


Appendix 
,»7T) — 2, VY) 
A(q/A)/h 
Ox, r) — Wx, 0) . oe 
en , dimensionless 
A(q/A)/h 
At(z, 7) — At(z,0) .. , 
—n , dimensionless 
A(q/A)/h 


Kr, dimensionless 


, dimensionless 


transformed function 


see equation (20) 
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ambient temperature. The transients are introduced by (a) 
suddenly increasing the plate temperature and (b) suddenly pro- 
ducing a uniform heat flux at the plate surface. The method of 
characteristics is employed to obtain the results which are re- 
ported as the approximate time required for the plate to reach 
steady state subsequent to the introduction of the transient. 
Anantanarayanan and Ramachandran [22] have studied the effects 
of traverse vibrations (75 to 120 cps) in a wire on the heat 
transfer to a turbulent air stream flowing parallel with the axis of 
the unexcited wire. Their results indicate that an increase takes 
place in the heat-transfer coefficient with an increase in both 
amplitude and frequency, and is correlated as essentially the ratio 
of vibrating to stationary Nusselt numbers versus mean vibra- 
tional velocity of air-stream velocity. 


Summary of the Present Problem 

The solution of the time response of the heat-transfer surface 
temperature and the fluid-surface temperature difference, at any 
point in a channel of this type of heat exchanger in which the rate 
of internal heat generation is an arbitrary function of time, results 
from a direct mathematical attack on the differential equations 
derived from the appropriate natural laws with the suitable 
boundary and initial conditions. The mathematical details are 
given in Appendix 2 and involve the use of the Laplace transform. 
While the solution is not algebraically simple, numerical results 
should be reasonably easy to obtain as all the functions necessary 
The 
authors have not rigorously proved the mathematical convergence 
of the complete solution to the steady-state values 


for such a computation are presented in graphical form. 


It appears 
from their own calculations in comparison with their experi- 
mental data, however, that such convergence is obtained. Asa 
special case, the solution also is given for an exchanger in which 
the heat generation is a linear function of time. These results 
are computed for a particular exchanger and compared with the 
experimental and theoretical results of references [1, 2] in which 
the transient was caused by a step change in internal heat genera- 
tion. Reasonable comparison is obtained. By way of this com- 
parison, it is shown how the case of linear time rate of change of 
heat generation may be used to compute the general case in which 
the heat generation is given by an arbitrary function of time, 
providing this function may be approximated piecewise by a series 
of linear variations. 


Analysis 


The physical system analyzed is shown in Fig. 1. This con- 
sists of a circular tube through which a coolant flows steadily, 
and in the solid walls of which energy is generated. In the steady 
state, all of this energy appears as a flow of heat at the interface 
between coolant and solid causing the coolant to increase in en- 
thalpy The out- 
During the transient, how- 
ever, both the tube-wall material and the coolant experience local 
temperature excursions, 


(and temperature) as it flows through the tube. 
side suriace of the tube is adiabatic. 


the magnitude of the former and their 
difference being evaluated here as a function of distance and 
time A circular tube has been selected for convenience in com- 
W here 


geometry is involved in the analysis, suc h 


parison of theory and experiment. a nonsymmetrical 
as a triangular flow 
cross section, the results given may be reduced to an area-to- 
volume ratio and, to a degree of approximation, can be used for 
other geometry for which the restraints on the solution apply.* 
The solution will be valid as presented for flow inside circular tubes 


and between parallel plates. 


’ The authors are indebted to Prof. T. F. Irving, Jr., and J. P. 
Hartnett of the Heat Transfer Laboratory, University of Minnesota, 
for pointing out the importance of this distinction 
Section [1, 2}). The nature of the approxin 
authors’ closure of reference [1 ]. 


see Discussion 
ation is discussed in the 
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Fig. 1 Simulated coolant channel 


The following assumptions are imposed on the solution: 

(a) The fluid temperature and velocity are constant across 
the flow cross section. 

(b) The duct wall temperature does not depend on radius. 
Valid for thin-walled metal systems. 

(c) Axial heat conduction is negligible in both fluid and duct 
walls and heat flows only to the coolant. This is a reasonable as- 
sumption when the Peclet number exceeds 100 [23]. 
surface is adiabatic. 


The outer 


(d) The heat-transfer coefficient is constant with length and 
time. 

(e) The fluid is incompressible and all fluid properties are con- 
stant. The flow channel has constant area. 

(f) Energy (heat) generation within the walls is constant with 
length but is time-dependent. 

(g) The temperature of the fluid entering the duct is constant 
and equal to & 


With these assumptions, the application of the First Law of 
Thermodynamics and the law of conservation of mass for an in- 
compressible fluid to the system shown in Fig. 1(b) produces the 
following two differential equations, one for the duct wall, the 
other for the fluid. 
Appendix 1. 


Details of this derivation are outlined in 


06 
OT 


Wall (0 —t) + ap,"(r) = a 


“— ol 
Fluid (@ t) = a; 


OT 


For a circular tube and in generalized terms, the constants are 
written: 


- pe,r; 


2h 


(3a) 


PCyT¢ 
a = u 


" (3t 
2h K ”) 


and p.(T) = Pz.” + P(r) (3c) 


where p,,” is the initial uniform volumetric heat-generation rate 
and (7) is the arbitrary time-variant volumetric heat-generation- 
rate function, having a zero value at zero time, which introduces 


the transient into the system. As shown in Appendix 2, both 
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functions Ox. r 1) and (2) 


two components, one being a steady-state part, 4(2, 0 and 6,(z, 


and ¢(2, r) of equations are split into 


()), the other of the transient contribution, 4(2, 7) and 6.(z, 7 

For purposes of transient analysis, it is necessary only to carry 
functions It 
ind #7, 0 & xr, 0). Hence 


the governing differential equations become 


out detailed mathematical analvsis of the latter 
will be noted that ¢,(2, 0 tr, U0 


Wall 


and (5) are operated on by the I iplace trans- 


wall 


form technique Phe 


transiormed equations in transient 
temperature (#@) and transient fluid-wall-temperature difference 
(At 0 t ure 


boundary 


integrated with the appropriate initial and 
conditions, the inverse transformation performed, and 
then recombined with the steady functions 4(2, 0) and 
4(2, 0 
mathematical attack on this problem produces solutions for both 


and At B. T 


physi al time 


state 
18 Outlined in detail in Appendix 2 The nature of the 
Hx, 7 in two domains of physical time: Case 1, a 
tlowing the introduction of the transient which is 
ind Case 2, a 
It will be noted that 
is that required for fluid particles flowing 


greater than zero but equal to or less than 2/1 
physi ul time equal to or greater than 2/u 
the physical time 
to traverse the distance z 


at a velocit The results for Cases 1 


and 2 for arbitrary timewise-variant internal volumetric heat 


generation oo T ire as follows 


Case 1 


a) W ill temperature Oz 


ind At L, 


0) represent the transient increase of wall te mperature 


In equations o and (7 
T Al 


and fluid-wall temperature difference at the location z 


the quantities 4 r.T 4 0 


r and 
The 


ire independent of the magnitude of #2, 0) or At(r, 0 


constants all are defined later 
Case 2 


Wall temperature 0 
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Al eT 


The quantities and functions given in equations (6 


(d) Fluid-wall-temperature difference Af(z, 7) 


At(z, 0) 


Ke 7 [ Af — sve ; 2k ol 
/0 4 Vu 


through (9 


are defined 


pe F J ws K 
pe,V K,, 


VW 


total heat capacity of wall material 


total heat capacity of fluid 


fluid time constant 


wall time constant 


(rr) = arbitrary time-variant volumetric heat-generation 


rate 


, and (9 


for the dynamic response of the wall temperature @(z, r 


Equations (6), (7 ), (3 represent the general solutions 
ind the 
initial 

heat 
Owing to their definition it will be noted that 


At(z, 7), 


and (7) or (8 ind (Y 


fluid-wall temperature difference Af(zx, 7) above thei 


values for an arbitrary transient in volumetric internal 
generation P(r 
the fluid temperature response U(r, T) = Oz, TF 
which may 
Should &(7 


change in volumetric heat source, that is, A(q/A)(4/J 


be formed from equations (6 


In equations (6 through (9) be set equal to @ step 
, it will 


be observed that the results of reference [2] are obtained 


Special Case of Transient in (-) Being Linear in Time 


\ simple but important case of the general solutions equations 


6) through (9), is that in which ®(7r) varies linearly with time 


This may be represented by the function 


T 
= A(q/A)A/) ( ), 
> 


are the volumetric heat-generation-rate 


11b) 


where Ap,” and A(q/A 
increment and the heat-flux increment, respectively, correspond- 
ing to any time T subsequent to the introduction of the linear 
transient 

For this case, as shown in \ppendix 2 the equations (6 through 
9) may be written in dimensionless form, as follows 


Case 1 0 < ru/z <1 


Wall temperature (x, 7) 
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Case 2 ru/z 


(g) Wall temperature #(z, 7 


~— ( 
MU +1 


15 


V 


In equations (14) and (15), q* and y* are to be interpreted as 
that value of q and y for which 7 corresponds to r* T z/u 
It will be noted Kr/M; Kr*/M. The 
functions necessary to the use of equations 14) and (15 ; Ws P. 
7"), Vols, 7"), Yo**(s, @*), and Wels, g*) are defined in Appendix 2, 
and have been calculated and are presented in graphical form in 
vs, a” 


y*), and ¥(s, g*) have been introduced in references [1] and 


that q hence q* 


3, 4, and 5 ) is a new function; yY.**(s, q*), 


Discussion of Theoretical Results 


Owing to the relative complexity of equations (14) and (15), it 
is simpler to examine equations (12) and (13) for the effect of the 
various physical parameters on the response of the wall tempera- 
ture, Wz, r 


At(x, r 


The constant M has been defined as 


6&2, 0) and the fl 1id-wall-temper ature difference 
A(z, 0) for the case of a linear transient in p,” 


Vv a, ei (Ply). (r,? (pe,\ - 


2 


ay (pe,) r3 pc,V 


total heat capacity of wall material 


total heat capacity of fluid 


Hence M is a heat-capacity ratio, dependent upon geometry as 
well as upon materials and is the quantity which introduces these 


properties into the equations. For the usual materials and 
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geometry it may be expected that M will vary from near zero for 
water-steel systems with thin metal walls to something in excess 
of 1000 for similar gas-steel systems. 

The effect of the heat-capacity-ratio parameter M on the wall- 
temperature response (for the first time domain) may be seen by 
an examination of equation (12) which may be regarded as a di- 
For small values of M, 
equation (12) approaches [1 + (1/2)Kr]r/r7o while at larger valuee 
of M, this quantity rapidly approaches zero. This indicates that, 
at least in the first domain (0 < ru/z < 1), the magnitude of the 
wall-temperature response is greatly diminished for a fixed Kr, 


mensionless wall-temperature response. 


for systems having a large value of M as compared with those of 
smaller value of M 
nature of the imposed transient. 


This result is a direct consequence of the 
That is, with the transient re- 
sulting from changes in the rate of internal generation of energy 
within the wall, the magnitude of the increase in wall temperature, 
at a fixed Kr and A(q/A)/h, will be inversely proportional 
to heat capacity (pc,V),, of the wall. Henceas M is increased the 
heat capacity of the wall relative to the fluid is increased and the 
expected reduction in the response of @ results. 

A similar result may be obtained for the response of the fluid- 
wall-temperature difference, also in the first time domain, by an 
examination of equation (13). For systems with small magni- 
tudes of M, that is, heat capacity of the wall very much smaller 
than that of the fluid, equation (13) approaches r/7». Hence for 
this case, Af(z, 7) — At(z, 0) approaches [A(q/A)/h]r/7o. On 
the other hand, for systems having very large M, equation (13 
approaches zero and Af(z, r) approaches Af(z, 0). 
approach to a condition of complete lack of transient which, like 


This is an 


that already mentioned for 6(z, Tr), may be expected from physi- 
cal considerations. 

The asymptotic limits of the temperature-response functions 
for values of M equal to zero and infinity are given for the first 
time domain in Table 1. 

Asymptotic limits for the second time domain have not been 
evaluated for the case of M = 0 owing to an indeterminate form 
of 0/0 in g* when r = z/uand r* = 0. Equations (14) and (15) 
indicate that each of the functions 8, 6, and @ has a limit of 0 for 
infinitely large values of M, as would be expected. 

As is shown in Appendix 2, the functions y(r, y*), Y2**(s, ¢*), 
Ys, g*), and Yols, g*),/which are required in equations (14) and 
(15) to describe the heat-exchanger response, are presented in 


9 
Figs. 2, 


“* 


3, 4, and 5 as various functions of the dimensionless quan- 
tities r, y, s, and q. 
written in the following form 


wens (7), (F) > ¥(3) 


It will be noted that these quantities may be 


h ( i , 
= )r = K, 
(pe, | 


u 


Table 1 Asymptotic limits of temperature-response functions for 
values of M = O and M 


First time domain 0 < 


(143 


Function 


O(x, r 6(z, VU 
A(q A) h 
Al(z, r) — Al(z, O 
A q A) h 
t(z, tr) — Kz, 0 
A(q/A)/h 
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Hence in the second time domain, ru/z > 1, the response of 
the heat exchanger also is governed by the product of a time con- 
stant for the fluid K or the wall K,, and either a physical time 7 or 
a “lag time’’ z/u. 
the area A are common to both wall and fluid, it may also be 
shown that A/K, = M and that sq = yr 


Because the heat-transfer coefficient h and 


Method of Superposition of Linear Transients 


The results for the linear transient given in equations (12), 
(13), (14), and (15) may be used to compute the response of 8, ¢, 
and At resulting from any other type of transient in heat-genera- 
tion rate providing this transient can satisfactorily be represented 
by a series of piecewise linear combinations. This technique 
follows from the method of superposition and is possible owing to 
the linearity of equations (47) and (48), which may be written 

6. function of Fis, Fis, Fis } = v . 
aybo/b,2 | from equation (47 | p> _— u 


t 


Ate function of F's, F 
aybe/b;? from equation (48 


where 


a p) = Di(p) + Pol P) 


eae 
> (Pp) (20 


¢,(p) = transformed form of an energy-generation transient 


having a linear form 


If one were to solve the equations (18) and (19) separately for 
each ¢,(p) and then superpose the resulting solutions, one ob- 
tains the solution to equations (47) and (48) in which the transient 
in energy generation ®(p) is made up of a linear combination of 
linear functions ¢,(p), in accordance with equation (20) 

\s an example, consider the transient shown in the inset of 
Fig. 6, labeled curve (c). In this case the energy transient con- 
sists of a linear increase until two seconds (7)) have elapsed, after 
Such a 
transient may be represented by an energy-generation rate from 
0 to 2 see by 


which the energy-generation rate remains constant. 


oir) = ap," (£) 


l 
T p* ) 


another energy 


and 


At 7 equal 7» (2 seconds generation rate function 


is introduced which is written as 


Ol T = 


—Ap,,” ( 


l l 
ww ie 
To p 


Then subsequent to rT = 7 


and o{p) = (23) 


). PT, 


the resultant condition of energy 


generation, which is constant and equal Ap,", is ex- 


¢ v4 
LO Pre 


- Ap." (=) Ap “(7 *) (24) 
T T 


= Ap.,” 


pressed as 


Pi P) = Gi(p) + Exp) = Ap.,” 
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t(x,2)-t(x,0) ,0(x,2)-@(x,0) °F 





160}— 


V*2.76 ft/sec 
4q*9.78 Kw 
tor 61° F 


Momentum Analogy Correction (1,2) 
* -_ 
Heat Balance 
Temperature (i,2)| 
——— Theoretical Responses 
Measured Responses (/,2) 
\UStep-Change in q) 


t(x,r)-t(2,0) ~ 





3 4 
t SECONDS 


Fig. 6 


Since e~ ?’s does not change the solution but shifts it byan amount 
To, it would be sufficient to solve the problem for the case (1/7o) 
(1/p*) and plot the resulting temperature curves for @, t, and At. 
Then after r = 7) subtract the same curve (starting at rT = To) 
The 


initial conditions are identical for both cases because the second 


from the previous one which begins at the origin (r = 0). 


solution (¢2) also actually starts from the origin (7,= 0) but in the 
region 0 < rT < 7) it is exactly equal to zero. 

The case of a linear increase in heat-generation rate, shown as 
curve (c) in the inset of Fig. 6, was computed and compared with 
both experimental data and theoretical results corresponding to a 
step increase in heat-generation rate taken from references [1, 2]. 
The final heat-generation rate and the initial conditions in the 
computed result were assumed to be the same as those correspond- 
ing to the experimental data. The results of the computation 
are shown in Fig. 6 for both @(z, rT) — @(z, 0) and At(z, r) — 
At(z, 0), and correspond to the transients shown as curves (a), 
(b), and (c) in the inset. Reasonable agreement is obtained, the 
computed response for case (c) lagging the step response, as would 
be anticipated. Also the asymptotic approach of ¢ to its steady- 
state value is identically the same for both transients and the 
experimental data, also as would be expected. 
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APPENDIX 1 
Derivation of the Basic Differential Equations 


The system analyzed is shown in Fig. 1 and consists of a con- 
fluid 
is pumped steadily and in the walls of which energy is dissipated. 

With the 
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stant-diameter circular tube through which a 


(coolant) 


assumptions and conditions as indicated in the 


NOVEMBER 1959 


analysis section of the paper, two differential equations, one for 
the wall and one for the fluid, may be derived from the First Law 
of Thermodynamics and the law of conservation of mass for an 
incompressible fluid [1 


re) | ph u 


where 


= flow wetted cross-sectional area (mr,? for circular tube) 
= wetted perimeter (27r, for circular tube) 


= cross-sectional area of duct wall (4,r? — mr,? for circu- 
lar tube) 
Px,” = initial volumetric heat-generation rate 
#(r) = arbitrary transient in volumetric heat-generation rate 
em 


internal energy per Ibu in absence of the effects of 
motion, capillarity, and so forth 


Equations (26) and (27) may be reduced to the following 


pea ( ot ) 
he . Or 
pe,T; ( ol ) 
Dh ; or 
( of r ( ol 
K Or K Ox 


is for circular-tube geometry and (28c) is for 
generalized geometry. 


(28a) 


(28c) 


l-quation (28b 
Both are equivalent to equation (2) in the 
body of the paper. 

For the wall, equation (27) may be rearranged to 


’ ee \ RY) 
. is...” + Pir | = : (pe, > ( ) (20a) 
he . he or 


— 1 (20 
4 = on” 4 P(7r)|) = — (29c) 
(pce, )wK » K, \ or 


Equation (296) is for circular-tube geometry and (29c 


generalized geometry. 


is for 
Both are equivalent to equation (1) in 
the body of the paper 

At this point, it is convenient to partition both temperatures 
Wx, r) and (xz, T) into transient and steady-state components. 
Thus we define 


Ox, r) = A(z, 0) + 0.2, r (30a) 


tz, T) = (2, 0) + 4(2, Tr) (306) 
where 6, and ¢; are the steady-state and @, and &, the transient 
components. For the steady-state condition, 6; and ¢, satisfy 
the following differential equations which are equations (28c) 


and (29c) with all time derivatives and P(r) set equal to zero, 
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(32) 


(33) 


6,(z, 0) = ty 4 - } - 34) 
where (q/A )» is the initial steady-state heat flux at r = 0. 
it will be noted that @(z, 0) = 6,(z, 0), 
“O, rT) = t,(0, Tr) 
and (2: 


Since 
(x, 0) = &(2, 0), and 
= t, we have the following conditions on 6 


Substitution of equations (30) into eq and (29¢ 


the result gives 


iS 

K Or 
l ( 002 
A OT ) 


the same as equations (4 


lations (28 , and 


and (32 


fen l (= ) 
= = or 


~ &r 
PCy wh, 


subtracting equations (31 from 


These equations are and (5) in the 


main body of the paper and are subje t to the initial and bound iry 
conditions given in equations (35 
The final solutions to this problem ure expressed is OC xr, 


Ox, 0), and tx, Tr t(x, 0) which are identical with @.(z, r 


to(r, T) Since, as was noted, #(z7,0) = @ an r, 0 


In this way, the solution (2, 7 A(x, 0 x, 7), for example, 


is the transient increase in @ above an arbitrary initial value 


(x, 0) and is independent of the value of 6 0 


, 8o long as the 


system conforms to the requirements ol 


APPENDIX 2 
The Integration of Equations (4) and (6) 


lquations 4) and (5) may be operated on by the Laplate 


operato! with respect to time T to produce the following 


transformed equations: 


el 


ry 
o 


The boundary 


to the results of Appendix 1 are: 


and initial conditions to be imposed according 
t.(z, 0) = 10a 
6(z2, 0) = 10h 


Hence equations 38) and (39) become 


at 
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Solving equations (41) and (42) for i,(z, p) gives 


dl, 
(43) 
dz 


The solution of this linear first-order differential equation, with 
the boundary condition 2,(0, p) = 0, is 


where 


tT G2p 


Substituting equation (44) into equations (41) and (42) 
solving for 62, gives 


Fy p oP P)F ep (47) 


In like manner, the fluid-wall-temperature difference 


Ate i 
is found to be 


Al 
ayb 


4,2 


where b, as + a; and b 2a 


the functions Fy, F ‘, Fy, and Fi. 
are Laplace-transformed functions in the variable p which have 
corresponding original functions Gi, Guy, G 


In equations (47) and (48 
, and Gi, in the varia- 
ble 7 obtained by performing an inverse transformation on the 
functions F. This operation is usually outlined in most text- 
books on operational mathematics, such as Churchill [24]. Ex- 
tensive tables of the functions F(p) and G(r) may be found in the 
and Erdélyi [26]. The La 
place-transformed functions F and the inverse-transformed func 
tions G appropriate to equations (47) and (48) are 


works of Campbell and Foster [25] 


listed in the 
following table: 
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) 


may be 


The product 
treated by the 


of transformed functions, such as F yi 
method of convolution, a technique of the inverse 
transformation for product of functions 


Churchill 24 


This is discussed by 


chap. 2, p. 36 Hence the transformed equa- 


tions in the physical domain of z and rT, which are the solutions to 


the differential equations (4) and (5), and for boundary and 
initial conditions of equation (40) are found to be 


Case 1 0 « 
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Restriction 


(7) = arbitrary time-variant volumetric heat-generation rate 
Equations (49), (50), (51), 


and (52) are the same as equations 
6), (4 8 


, and (9) in the main body of the paper 


Special Case of Transient in #(r) Being Linear in Time 


In the event P(r) be represented by a linear function of time, 
such as Ap,"(r/r (50 51 


\ ’ 


, equations (49), , and (52) become 


Case 1 


_ 59 
a 

The functions in these equations are defined. The function 
y.**(z, T) is identical with two functions which have been com- 
puted and presented by Rizika [11]. 
duced here in equivalent form in Fig. 4. The remaining functions 
¥(z, rT) and y¥,(z, 7) have been presented in reference [2] 
here in equivalent form in Figs. 2 and 3 
with this paper and is given in Fig 5. 


The functions are repro- 


given 
¥(z, 7) is a new function 
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DISCUSSION 
G. M. Dusinberre* 


While this paper is admirable 


ithematics, it is opel 


problem 
2 The 
that heat generation is constan 


many applications to nuclear r 


m is readily solved 


ithors’ reference [12], with the addition 


4 Department f Mechanica Eng 


State University niversity Park, Pa 
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term. In fact, a simple hand calculation will easily reproduce 
the result shown in Fig. 6 of the paper. 

Table 2 of this discussion shows data from Fig. 6, from the 
authors’ reference [1], and certain data which were not stated 
by the authors and had to be assumed. (The quantitative agree- 
ment of the results depends on a lucky guess as to the specific 
heat of the material. Qualitative agreement would have béen 
attained in any case.) 


Table 2 
From Fig. 6 
h = 
lmean = 


From reference [1 
Length = 35 in. 
Material—stainless steel 

Assumed 
Cp = 0.120 Btu/lb deg F 


Table 3 
Flow area, Ay = 1.364 XK 1078 sq ft 
Heat-transfer area, A, = 0.382 sq it 
Surface coefficient, A = 700 Btu/hr sq ft deg F 
Flow rate, w = 2.76 X 3600 X 62.2 X 1.364 K 10-3 = 844 lb 
per hr 


Capacity rate, V = 844 X 1.0 = 844 Btu/hr deg F 

Heat generation, q = 9.78 X 34138 = 33,400 Btu/hr 

Steady-state temperature rise = 33,400/844 = 39.6 deg F 

Weight of tube = 0.526 Ib 

Heat capacity of tube, W = 0.526 X 0.120 = 0.0631 Btu/deg F 

Weight of water in tube = 0.248 |b 

Heat capacity of water, C = 0.248 Btu/deg F 

Surface conductance of tube, Al = 0.382 Kk 700 267 Btu/hr 
deg | 


Table 3 of this discussion shows certain preliminary calculations 
of a routine nature. Table 4 shows the actual calculation of the 
first 2 sec of transient (a) of Fig. 6. If the results of column 4 
and column D are plotted on that figure, it will be seen that the 
discrepancy is negligible. 

The preliminary calculations took about an hour. The time 
required for the calculations in Table 4, by slide rule, was 24 
min. This seems to substantiate the claim that the authors’ 
method is too difficult a way to solve a relatively simple problem 

It must be noted that an extremely crude subdivision has been 

d, and in spite of this the results are comparable with those of 
the paper. The numerical analysis may be refined to any degree 
which is warranted by the precision ol the data At the same 
time the computer programming is elementary. All the detailed 
calculations are explicit. No evaluation of series is required. 

It was observed that condition (f) is quite restrictive in nuclear- 

ctor applications. But in the numerical analysis there is no 
lifficulty in introducing information as to the actual distribution 
of the heat generation. 

To demonstrate this, Table 5, herewith, shows the result of 
calculations for the same total generation as the foregoing, but 
vith 45.6 per cent ol the generation in the middle third of the 
tube. The calculation runs for 3 s« c, and the steady-state data are 

Iso shown Column D is not computed, as the highest wall 


ture is found near the center. 


For a crude solution, assume three subdivisions, as in Fig. 1 of 
reference [12]. 
H 267 /3 89 W 0.0631 /3 0.0210 
( 0.248/3 = 0.0827 q 33,400/3 11,130 
H /2\ 89/1688 < ] 
Ar W/H 0.0201 /89 0.000275 hr = 0.99 sec 
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Point A 

Coefficients 0.196 
0.608 
0.196 
49 


0.000008 | 


LO.SOO hye 
I, 2, 3, 4¢ fluid 
\ B, C desig ale teri linte points the 


1) designates ou 


0.100 
LO SOO 


O.805 


0.106 
10.800 


Org U.0O0OS 


O.S805/ 0.100 + O.O05/ 


O.100/ O.1O0/ tr O.OOSt, + 49.0. ete. 


Use linear extr ipolation tor ¢, 
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129 


Table 5 


11,130 
19.0 deg 
0.0210 & 10.800 


Authors’ Closure 

We ippreciate the contribution Profe SSor Dusinbe rre 
in his discussion of our work but regret that he has not 
this opportunity to show in a broader sense the re! 
umong various computational procedures lor engineer! 
lems Our problem and analytieal s lution and Professor 
Dusinberre’s method provide an excellent vehicle for doit g this 
It should be pointed out that Professor Dusinberre has lirected 
his discussion entirely to the results of our two previous papers 

1,2). The subject of this present paper, which is an extensior 
of our work reported in 1, 2] ind signifie intly broadens it, Is no 
specifically discussed. 

It is possible that Professor Dusinberre has overlooked the 
principal purpose of our work. It was not our intent to determin 
the dynamic response of a '/.-in. ID stainless-steel tube hav ing 

-in. wall thickness and a length of about 3 ft, as Professor 
Dusinberre appears to suggest, nor to solve any othe! specific 
heat-ex¢ hanger problem. father, it was to obtain a rigorously 
complete and resolutely exact analytical solution to i model 
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which closely resembles the physical system corresponding to an 
the heat ex- 
We know of no such pre- 


important engineering problem of current interest 
changer having internal heat sources. 
vious published solution, although the literature is quite large on 
problems of a similar nature. Mathematical solutions continuous 
in space and time have the advantage of exactness and are subject 
to parameterization, as is well known 

We shall not attempt to defend our solution “from the point of 
view of practical engineering’ as this type of criticism is pretty 
much a matter of experience, the demands of one’s professional 


work and his personal taste. What is practical to one person may 


be impractical to another, ete However, within our own expe- 
rience we have found the use of applied mathematics to be of great 
But this does not 


constitute our defense of applied mathematics as it requires no 


value in the solution of engineering problems 
defense, at least not by us. We wish to be careful, however, to 
emphasize that we also have found numerical methods to be power- 
ful tools in engineering practice and we have used them in in- 
stances where a solution to a problem could be obtained in no 


other way. We 
on the importance of numerical methods and we fully support his 


are in full agreement with Professor Dusinberre 
€ mphasis of this In fact, as we have ittempted to make clear in 
the introduction to this paper and in both the previous papers 
1, 2], we recognize that a numerical solution to a problem of this 
sort is often convenient and practical 
had to 

dures which required the use of an IBM 650 « omputer, as was indi- 
cated Furthermore is we 
Dietrich [15 


Actually, in obtaining our 


‘analytical’’ solution we resort to numerical 


proce- 


stated 


already have Silvers and 
solution of 


Hence, 


this, as Professor 


have used numerical methods in the 
the response of a specihec heterogeneous nuclear reactor 
it does not seem necessary to redemonstrate 
Dusinberre has done, although it probably is of some value to re- 
numerical methods. We 


thought we had freely recognized that in the paper 


emphasize the use and importance of 


Like analytical methods, numerical procedures require ho 
defense nor apology Both are well established, recognized tech- 
niques used often in concert to bring about a solution to an engi- 
neering problem. We would champion the knowledge of both 
for engineers, not one to the exclusion of the other. 

Professor Dusinberre states that our solution is much too diffi- 
eult 


found and the necessary functions given it is irrelevant to be con- 


The task of obtaining it was great but now that it has been 
cerned with the magnitude of our task, since to use our results 
does not require one to repeat the efforts we put into the problem 
We al- 


ir results is quite simple, 


Certainly one does not have to evaluate 
ready The use of o 
requiring a knowledge only of algebra 


a single series. 
have done that 


while considerably more 


knowledge and experience is required to solve the same problem 
by numerical methods 

In a sense Professor Dusinberre’s criticisms are not reasonable 
He states that by using even a crude subdivision his numerical 
computations provide a result essentially equivalent to ours and 
imphies, therefore, the futility of all but a numerical solution to 
this problem. But this is hardly a proper evaluation of this situa- 
tion. It isimportant to recognize that in forming this conclusion 
he relies on our exact solution as a criterion for the completeness 


That 


portance of the essential agreement between his solution and ours. 


of his own numerical solution is why he stresses the im- 
However, the question now arises as to what one would do if an 
Would one be 


network as crude us dividing a channel length into thirds? We 


exact solution were not available. satisfied with a 


doubt it, if a solution is being sought to a significant engineering 
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problem. Does Professor Dusinberre’s result now mean that it is 
unnecessary to divide a channel into any subdivisions smaller 
than one third? Probably not. 
our exact solution that he is able in a single attempt with a crude 
Without 


the results of an exact solution available considerably more effort 


Hence, it is only because he has 
network to obtain a solution with which he is satisfied. 


would be required to obtain the solution by numerical means than 
the one hour and 24 minutes it took Professor Dusinberre. It 
would be our guess that something like two or three more finer 
networks would be tried to establish the convergence of the solu- 
tion to its true value. In such a case the time required would be 
In this 


regard it is important to remember that Professor Dusinberre’s 


considerably greater even for an experienced worker. 


solution has been produced by one of the highest authorities on the 
subject of numerical methods and his one hour and 24 minute 
While 


this is somewhat speculative, it would seem to us that an engineer 


computation time cannot be taken as par for the course. 


with average experience in these matters could well take some- 
thing more like two days to work out a satisfactory numerical 
solution to this type of problem and still give due consideration to 
the convergence and stability of his result 

We believe that what this discussion really has shown is the 
interplay between an exact analytical solution and a numeri- 
cal solution to the same problem. The former often is of great 
value in computational work as a check solution on a numerical 
program, acting as an important point of reference 

It is possible to bring about a simplification of our results in 
another way. 


Taking equation (48) as an example, we find that 


the complex form of the final solution arises from the term 


However, by the use of reduction techniques given in reference 
5] this term may be approximated as follows: 


If this approximate form is substituted into equation ( 18) one 


finds: 


Ap,” ; : 
transformation of the 


l 
- then an inverse 
To p* 


loregoing gives the following two results 


If we let ca 
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Aer ae 


rhis approximation, which also could be applied to equation (47 
for 6, and reduced to @,(z, 7), requires the calculation of three z- 
dependent parameters and the use of three exponential functions 
However, only the use of algebra is nece ssary and the time re- 
quire d to obtain a solution would be no k nger than that for Pro- 
fessor Dusinberre’s method. 
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Heating and Cooling Air and Carbon Dioxide in 
the Thermal Entrance Region of a Circular Duct 


HELMUT WOLF 


Senior Research Engineer, 

Research Center, Rocketdyne, 

A Division of North American Aviation, Inc., 
Canoga Park, Calif. Assoc. Mem. ASME 


With Large Gas to Wall Temperature Differences 


The local heat-transfer characteristics for air and carbon dioxide in the thermal entrance 
region have been determined experimentally and are compared with predicted values 
computed from an extension of the theoretical analysis due to Deissler. 


The expert 


ments were conducted in smooth round tubes having a fully developed turbulent velocity 


profile and a uniform temperature distribution at the entrance. 


The boundary condt- 


tions of uniform heat flux for heating the gas, and wall temperature constant and varia- 
ble with axial distance for cooling the gas, were investigated with large temperature dif- 


ferences between the gas and the tube wall 


For heating, the experimental results 


yielded one per cent thermal entrance lengths ranging from 11 to 27 diameters over the 


bulk Reynolds number range from 50,000 to 246,000 


For cooling, one per cent thermal 


entrance lengths ranging from 12 to 26 diameters were obtained for bulk Reynolds num 


bers ranging from 17,000 to 218,000. 


The agreement between theory and experiment 


was favorable in most cases. 


Introduction 


= A FLUID FLOWING ADIABATICALLY in a circular 
duct enters a region having a wall temperature different from that 
of the fluid, the temperature difference may be 
matically 


treated mathe- 
as a discontinuity. At that cross section, termed the 
entrance, the temperature gradient in the fluid layer near the wall 
The resulting local heat-transfer coefficient 
and Nusselt number are in theory likewis« 
rapidly in the direction of fluid flow 


is, in theory, infinite. 
infinite and decrease 
In the design of heat- 
transfer apparatus, it is of considerable interest that quantitative 
data be available regarding the rate of heat transfer to or from 
the fluid as a function of the distance from the entrance. 

This paper summarizes the results of an analytical and experi- 
mental investigation of the heat-transfer characteristics of air and 
carbon dioxide in turbulent flow in the entrance region of a 
smooth round duct with large temperature differences between 
the gas and the wall. In the analytical investigation the fluid 


was assumed to have a fully developed velocity profile and a uni- 


! The research reported in this paper was performed at the Jet Pro- 
pulsion Center, Purdue University, and was sponsored by Project 
SQUID. 

Contributed by the Heat Transfer Division and presented at the 
Semi-Annual Meeting, St. Louis, Mo., June 14-18, 1959, of Tue 
AMERICAN Society oF MECHANICAL ENGINEERS. 

Note: Statements and opinions advanced in papers are to be 
understood as individual expressions of their authors and not those 
of the Society. Manuscript received at ASME Headquarters, March 
4, 1959. Paper No. 59—SA-15. 


form temperature profile at the cross section where energy ex- 
change was initiated. In the experimental investigations, a duct 
length 30 diameters long was provided upstream to the test sec- 
tions to insure a fully developed velocity profile (with a uniform 
temperature profile) at the entrance [24].? 

Theoretical heat 
region are presented in references [1, 2, 3, 


transfer in the entrance 
4, 5, 6) for the flow of 
7, 8] for flow between flat 


calculations of 


fluids in tubes and in references [3, 
plates with the wall temperature constant with axial distance 
teferences [5, 9, 7, 20] discuss the calculation of heat transfer in 
the entrance region for the constant heat flux boundary condition, 
but assumed that the physical properties of the flowing fluid did 
not vary with 


temperature. Only the analyses presented in 


references [5, 9] permit including the temperature dependence of 
the fluid physical properties. Experimental determinations of 
heat transfer in the entrance region are reported in references {10, 
11, 12, 13] for Prandtl numbers ranging from 0.7 to 350 and for 
temperature differences between the fluid and the wall ranging 
from 8 to 390 R. Cooling experiments with large temperature 
differences are described in reference [14] at a Reynolds number 
of approximately 10,000; however, the variation of the local heat- 
transfer coefficient with axial distance was not determined. 
The investigation reported herein encompasses the following 
conditions: |, — 2,; > 400 R, Reynolds numbers ranging from 
17,000 to 246,000, Prandtl number 0.7, and for heating and cool- 


ing the gases. 


2 Numbers in brackets designate References at end of paper 





Nomenclature 


= inner surface area of the test h, = 
section, ft? 


> } - 
= cross-sectional area of the duct a 


ft2 efficient, Btu/hr ft? R 
= thermal conductivity, Btu/hr 


specific heat at contant pres- , 


sure, Btu/lb ftR 


= heated length of the tube, ft Q, 


inner tube diameter, ft 
heat balance, per cent 
gravitational constant, 32.17 a = 


{t/sec? 
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local heat-transfer coefficient, 
3tu/hr ft? R 
= asymptotic 


= one per cent entrance length, 
diameters 

five per cent entrance length, 
diameters 


Nu = 
Nre 

Np, Prandtl number, 3600gc,u/k 
Qai 


Nusselt number, hD/k 

teynolds number, (D/A .g)w,/pu 

heat-transfer co- 

energy per unit time passing 
through A,, Btu/sec 

, energy per unit time absorbed 

by or extracted from the gas 

in passing through the test 

section, Btu/sec 


(Continued on next page) 
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Fig. 1 Schematic flow diagram of heating system 


Experimental Apparatus heating and cooling experiments was preceded by a 30-diameter 
Two different sets of apparatus were employed for the heating isothermal starting length and followed by a 10-diameter isother- 
and the cooling experiments. Each apparatus was designed so mal length. The surface temperature of each test section was 
that either air or carbon dioxide could be admitted to the test eec- measured at 16 stations along its length by chromel-alumel ther- 
tions without interrupting the gas flow. The gas-flow rate was mocouples capacitance welded to the surface. The temperature 
measured with calibrated sharp-edged orifices in conjunction with of the gas in the inlet and exit mixing chambers was measured 
appropriate water manometers. Simultaneous control of a dome with chromel-alumel thermocouples carefully calibrated against 
loaded pressure regulator and a variable area orifice downstream * standard platinum versus platinum-10-per-cent-rhodium ther- 
from the test section permitted varying the gas-flow rate in the mocouple whose emf-temperature relationship was determined by 
test section while maintaining the pressure at its entrance at a the National Bureau of Standards, Washington, D. C 
selected value Apparatus for Cooling Experiments. Fig. 2 presents the flow diagram 
Apparatus for Heating Experiments. Fig. | presents the flow dia- of the system employed in the cooling experiments. Suitable 
gram for the apparatus employed in the heating experiments baffles were provided in the mixing chambers for thoroughly 
Three different test sections, fabricated from Havnes Allov 25 mixing the gas prior to measuring its temperature with calibrated 
L605) tubing, were employed in the heating experiments. All of shielded platinum versus platinum-10-per-cent-rhodium ther- 
the tubes had an ID 0.569 in. and an OD = 0.629 in.; their mocouples. The test section was fabricated from type 347 stain- 
L/D ratios were 21.1, 40.3, and 59.75. Each test section for the — less-steel pipe, precision ground to 0.535 in. OD and internally 


Nomenclature 





energy per unit time and area y = co-ordinate perpendicular to Subscripts 


passing through the wall at the duct wall, ft b = pertaining to or evaluated at the 


bulk temperature 

the wall, UV Tr. Po/ (My / Pw ] pertaining to or evaluated at the 
Patho €47/ Ex initial temperature 
heat-transfer parameter for 


location z, Btu/ft? see ’ nondimensional distance from 
radial distance from the center 

line of the duct, ft 
temperature at the point (2, y ares 
temperature of the gas at 2 <0, ; pertaining to or evaluated at the 
wall temperature invariant 


constant with respect to y wall temperature 


with axial distance, 
temperature = at the location 


' pertaining to the thermal boundary 
y 5, eV 7. Prwo/ CpQT wl layer 


nondimensional temperature = Ph btm tac . 
; t/t,)/B arv laver Superscripts 
velocity in the z direction at the : eddy diffusivity for heat and * = pertaining to a nominal value 
point (2, y momentum, respectively RB = corrected to a nominal value of 
nondimensional velocity, = mass per unit volume, lb sec? Reynolds number and 6 
r Vr, Pp. it* - a bar over the symbol denotes an 
co-ordinate parallel to the wall, = dynamic viscosity, lb sec/ft* average quantity over the length 
ft (x = Oat the entrance tT = shearing stress, lb/ft? of the test section, L 


168 / NOVEMBER 1959 Transactions of the ASME 





PSIG / 
1200 PSIG 
CO» CYLINDERS (c) 
WATER 


FLOW ORIFICE 


(G) 2 


INL MIXING CHAMBER 
(1) 


PROPANE 7 ~ 
TEST SECTION 
PRESSURE DROP 

MANOMETER 


—— TEST SECTION (K) 


HOT WATER- 


7SF WATER 


G 
aTeRs (™) 


== 


xIT 
CHamBer ‘+? 


WATER SOFTENER 
LEOGEND 
HAND LOADER 
PRESSURE GAGE 
CHECK VALVE 
HAND VALVE 
DOME LOADER 
METERING VALVES 


TOTAL FLOW RATE THERMOCOUPLES 


WEIGH SYSTEM 


Fig. 2. Schematic flow diagram of cooling apparatus 


nd was fitted with 
Water 


with low values 


precision bored and polished to 0.384 in. 1D 


a cooling J icket containing 12 separate coour 


ig sections 
vas employed as the coolant in the experiments 
The flow rates of the 


the individual sections were measured with 12 


of the wall temperature cooling water to 
ilibrated rotame- 
ters Fig 2 also shows the svstem « mplo ed iors ipphy Ing a Uunl- 
the test 
At values of the wall t« mperature greater 


form rate of cooling water at constant temperature to 
section vooling ja ket 


than 710 R, 


\ detailed dese1 iption ol the components « 


the test section was cooled with ai 
my] loved in the heat- 
ing and cooling systems and of the various calibrations is given in 
relerence 16). 

Heat Balance. ‘The heat balance e, stated in per cent, is defined 
by 
a) a) 


«Ai ‘ 


For the electrically 


computed as the difference between the energy Q 


a e 


heated test sections iantitv G0 


day WAS 
gene! ited by 
the current passing through the wall of the test section and the 
energy loss Q; by conduction through the insulation surrounding 
s assumed 


The 


average of the absolute values of the heat balances obtained for 


it. For the cooling experiments, the quantity Q4; wa 


Ai 
eq tal to the sensible energy increase ), of the coolant 


the heating and cooling experiments with air and carbon dioxide 


are summarized in Table 1. 


Table | 


Heating experiments; air 

Heating experiments; carbon dioxide 2 
Cooling experiments; air 3 
Cooling experiments; carbon dioxide } 


Average absolute values of the heat balance 

05 pe 
16 pe 
1) pe 
15 pe 


Experimental Results 


Local Heat-Transfer Characteristics, Heating. Since the 


heating 


periments were conducted in three different test sections, it 
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not possible im the experiments to maintain identical flow condi- 


Also, the local values of the heat- 


the parametet! 8B is » measure of the in- 


tions in all three test sections 
transfer parameter B 
tensity of the heat transfer and may be considered approxim itely 
proportional to the wall-to-bulk temperature difference) were not 
constant along the length of a test section ( onsequently, to 
local Nusselt 
numbers were corrected to a nominal value for the bulk Reynolds 


* and p* 


correction is presented in detail in reference [16 


correlate the results for the three test sections, the 


The method of making the afore-mentioned 
The wall shear 


isured friction pres- 


number \ Re 


was evaluated from the experimentally mx 
sure drop based upon the average bulk velocity com- 
ind 7, 174, 
part Il, gives a detailed derivation and discussion of the par imeter 
B The corrected values of the local Nusselt number Nu, 8 lor 
air and carbon dioxice 


of the gas 


puted from p in the test section. Reference [16], p 


are prese nted in Fig. 3 as functions of the 
with Vpro* 
3 that N 1,8 decreases continuously along the length of 


ixial distance r/D and B* as parameters. It is seen 
from Fig 
the test section, due to the increase in the thermal conductivity of 
the gas with the bulk te mperature, The latte: figure demonstrates 
that as B* is decreased from 0.04 to 0.03 the value of Nu," de- 
‘ The afore- 
mentioned trend is in agreement with that predicted by theory 
[15, 16). 


ing the physical properties of the gas at the loeal bulk tempera- 


creases more rapidly at the higher \ luc s of Ny 


The Nusselt number Nu, #4 was computed by evaluat- 


ture and correlated with the parameter 8. Deissler [5] recom- 
mends the parameter @ = q,/c,p;u,t; for correlating constant heat- 
flux data with the physical properties evaluated at the initial 
temperature ¢ The parameter > is employed in the discussion 


‘- 
of the thermal entrance length. The « xperimental values of the 
local heat-transfer coefficients for the heating experiments are 
presented in reference [16]. 

Local Heat-Transfer Characteristics, Cooling. Fig. 4/1 
local Nusselt number Nu, for air [Fig. 4(6) 


constant wall temperature as a function of the dimensionless axial 


presents the 


for carbon dioxide | at 
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uniform initial temperature profile, 


that for Ol,,/Ox greater than zero, the local value of ¢ 
the thermal « 


of the gas changes little along the length of the test 


does not 


decrease appreciably. Consequently, onductivity 
For 


»t,,/Or less than zero, the bulk temperature of the gas decreases 


I ipidly 


section 


so does the thermal conduc- 


with axial distance, and 


tivity of the gas ky The experiments showed that the loe il heat- 
Ap- 


h that the 


transfer coefficient h, also decreased with axia 


distance 
parently, the rates of decrease for h, and k, were su 
value for Nu, remained constant 

The local heat-transfer coefficient h, is defined by the following 
equation: 
> 


I:xpressing (q,,) by Fourier’s law and substituting the result into 


equation (2) yields 


h, = k(dt/dy),,/(t 


w 
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Fig. 4 Variation of local Nusselt number with x/D, Reynolds number, and wall temperature gradient for air and carbon dioxide; 
cooling 


Analytie inve rations of heat transfer in circular ducts have teference 16) presents tables of = Nu, Nu, Yu Nes, p, 2, 
demonstrated that for laminar and turbulent flow [15, 17, I8]h, 4, and ¢, as a function of x/D for the heating experiments and 


decreases rapidly in the entrance region and then attains a sub- tables of h,, Nu,, g., Vres, B, t,, and t, as a function of 2/D for 
stantially constant value. Reference to equation (3) shows that the cooling experiments conducted with air and carbon dioxide 

, to remain constant outside the entrance region (0t/Oy),, Thermal Entrance Length. The thermal! entrance length is defined 

ary in the same manner as (t, or constant k. On as the distance from the entrance to the cross section where the 

the other nd, if the wall temperature is not constant with temperature profile of the gas has no straight line portion perpen- 
axial distan: hen in equation (3) varies with axial distance dicular to the tube axis A precise determination of the entrance 

Consequently, it cannot be expected that length, as defined, requires a calculation or measurement of the 

onstant value. thermal boundary layer. When that is not possible, the en- 

the local heat-transfer coefficient h, asafune- trance length is determined as the distance required for either 


tion of the axial distance x/D for air [Fig. 5(b) for carbon dioxide } Nu or h, to reach a value constant with axial distance Figs 
for constant wall te mperature Figs. 5(a and 6) demonstrate that through 5 and reference [19] demonstrate that, for fluids having 
ifter the first 12 to 26 diameters, depending on the Reynolds ‘physical properties that vary appreciably with temperature, the 
number, the value of the heat-transfer coefficient h, is prac- local Nusselt number is not a reliable criterion for establishing the 
tically constant. Figs. 5(c and d) present h, as a function of 2/D value of the thermal entrance length. Accordingly, the thermal 
for the cooling experiments conducted with air and carbon dioxide entrance lengths reported herein were inferred from the local 
for positive and negative values of the wall temperature gradient. _heat-transfer coefficient h,. Because it is diffeult to de termine 
Figs. 5(c and d) demonstrate that for variable wall temperature precisely the point at which the local heat-transfer coefficient 
h, continues to decrease by a small amount for the entire length when plotted as a function of z/D, reaches a value constant with 
of the test section; h, decreases more for carbon dioxide than for 2/D (denoted by h,), it is convenient to define a “5 per cent”’ and 


ur. Ghee per cent” entrance length Ln and Li respec tive ly The 
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Fig. 5 Variation of local heat-transfer coefficient with axial distance, Reynolds number, and wall temperature gradient for air 


and carbon dioxide; cooling 


quantit the distance from the entrance to the cross section 
1.05h,; similarly, Ly 


trance to the cross section where h, 


where A 


. is the distance from the en- 


101A, 

Fig 6 illustrates the variation of La is a funetion of the initial 
Reynolds number Vy, with @ as a parameter for heating air and 
carbon dioxide. The results presented in Fig. 6 demonstrate that 


the values of La for air are significantly larger than those for 


carbon dioxide at a given Ny 
variation of Ly, as a function of the bulk 
with pb iis iu 


carbon dioxide with the wall temperature ¢, 
distance The 


Fig. 7 presents the 


Rey nolis number par imeter tor cooling air and 
constant with axis! 
values of L,, obtained for earbon dioxide are ip- 


The 


computed 


than those for ai 
theoretically 


the agreement between theory and 


proximately 
solid 


values for 


8 to 10 per cent smaller 
Figs. 6 and 7 are the 
length 
experiment is conservative in all cases except that for heating car- 
Table 2 


ind Ly, for the different experiments conducted by 


lime Ss in 
entrance 
bon dioxice summarizes the range of values deter- 
mined for Ly 
the writer 
Over the range of Reynolds numbers presented in Figs. 6 and 7 
the values of Ly); can be represented, with an accuracy of +10 per 
cent, by the following empirical equation: 


Ly = CNe.' 
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Table 2 Range of experimental values for the thermal entrance length 


Gas L an V; \ 

Air, heating CHF 9-16 60,000 to 50,000 
200 , 000 158 000 

CO., heating CHF 100,000 to 46 COO 
300,000 000 
Air, cooling CWT 000 
52.000 
Air, cooling VWT P : ) ‘ 20 000 
04.000 

CO., cooling CWT f : 23 000 
218.000 
24,000 
127.000 


CO., cooling VWT = 12-28 


CHF denotes constant heat flux. 
CWT denotes wall temperature constant with axial distan 
VWT denotes variable wall temperature with axial distan 
Estimated values. 


Values of the constant C are given in Table 3 for heating 
cooling air and carbon dioxide 

Table } 
perimental determinations of entrance length with the results of 


Prandt] 


summarizes and compares the results of other ex- 


the subject investigation For fluids having numbers 
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Fig.6 Variation of thermal entrance length with initial Reynolds number 
and ¢ for air and carbon dioxide; heating, constant heat flux 


Table 3 Values of the « fant C in equation (4) 
Temperature for 
evaluating the 
properties 

Initial 

Bulk 

Initial 

Bulk 





Heating or 
Gas cooling 
Air 
Air 
CO, 
CO, 


Heating 
Cooling 
Heating 
Cooling 


between 0.7 and 8, Table 4 indicates entrance lengths from 8 to 35 
diameters for Reynolds numbers ranging from 4000 to 250,000. 

The phy sical properties employed in the correlation of the ex- 
perimental data were those recommended by Nuttal, Benedict, 
Morey, and Masi in reference [21]. The afore-mentioned proper- 
recommended by Hilsenrath and 


ties are substantially those 


Touloukian [22]. 


Theoretical Analysis 

Of the heat-transfer analyses published in the open literature, 
NASA permits a 
readily obtainable solution for the case where the physical proper- 


only the work of Deissler and co-workers at 


ties vary strongly with the temperature. The method presented 
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Fig. 7 Variation of thermal entrance length with bulk Reynolds number 
and § for air and carbon dioxide; cooling, fant wall ¢ 





ure 
by Deissler [9] assumed the specific heat and Prandtl number con- 
stant with temperature, identical variations of the viscosity and 
thermal conductivity with temperature for the constant heat-flux 
boundary condition only with properties variable. 

Accordingly, the analysis by Deissler [9] was extended by 
employing an approximation for the variation of viscosity, ther- 
mal conductivity, and specific heat with temperature that closely 
agreed with the actual variation for air and carbon dioxide, and 
by investigating the constant wall-temperature boundary condi- 
tion for the case of cooling the gases. 

The assumptions employed in calculating the velocity and tet 
perature profiles for fully developed turbulent flow employed in 
the analysis are given in references [15 and 16]. 

Boundary-Layer Equation. By making an energy balance on the 
boundary layer in the manner of Eckert [24], it is possible to re- 
late the boundary-layer thickness 6, to the axial distance from 
the point where either heating or cooling of the fluid is initiated. 


Thus 
d “s 
c,git — pu(r, y)dy dx 
dr /0 


du Vy dx = 


Table 4 Summary of experimental determinations of entrance lengths in tubes; uniform initial temperature 


profile, fully developed initial velocity profile 


Prandtl 
Boundary* number 
Investigator condition 
Boelter [10] CWT 
Aladyev [11] CWT 
Hartnett [12] THF 
JHF 
WT 
WT 
THE 


Chollette [13] 
Wolf 
(Purdue) 


0.7 
UHF denotes uniform heat flux 
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lw — bbl, 
degrees R 


48 


300-1200 
500-1000 


entrance 
teynolds number length 
N Reb (one per cent) 
27 ,000—53, 000 10-17 
4000-60 , OOO 35-20 
17 ,000-89 , 000 8-15 
5500-46, 600 29-9 
8000-12 , 500 10-15 
23 , 000-218 ,000 12-26 
50, 000-246 , 000 11-27 


140 
150 
8-72 
390 
140 


* CWT denotes wall temperature constant with axial distance. 
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Fig. sy Variation of local initial Nusselt number with axial distance, 4, and initial Reynolds number for air and carbon dioxide; 


constant heat flux (heating) 


For the case where the wall temperature is constant with axial 


distance qg, is a function of z. Integration and nondimensional- 


ization of equation (5) in terms of u* and yields 


be 
directly for the constant heat-flux boundary condition; in that 


Constant Heat-Flux Solution. Equation (5) may integrated 


case q, is not a function of z. Integrating equation (5), dividing 


both sides by 2r,?, and nondimensionalizing the result gives 


“ 
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To distinguish between the case where q,, is not a function of z, 
heat flux brackets 
Deissler [9] has proposed a parameter @ as a measure of the 


, constant z/D is not enclosed by 


1.e 


intensity of the heat transfer under the afore-mentioned condi- 
The parameter @ is defined as follows: Thu 


tion 


y 
Vu 


CUP iuit 


The parameter @ may be related to the parameter £6 | 
that for a constant area duct pju; = p,u, and t;* = t3* 
sult is 


noting 
The re- 


1 — Bt,* 
u,*(1 — Bts* 


ti) 


After solving equation (7) heat transfer and flow parameters were 
calculated for selected values of 8 as a function of z/D. The re- 
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Fig. 9 Variation of local bulk Nusselt number with axial distance, 3, and Reynolds number for air and carbon dioxide; constant wall tem- 


perature (cooling) 


A detailed 
description of the procedure is presented: by Wolf and Zucrow in 


sults were then cross-plotted for constant values of @ 


. 99 
relerence (Zo 


The derivation and method of calculation of the bulk velocity, 
bulk temperature, Nusselt numbers, and Reynolds numbers, and 


a omple te de scription ol the numerical methods employed for the 
solution, and of the program for the computer is given in ref- 
erences 16 ind 4} 


Analytical Results 


Values of the bulk Nusselt number Nu, were computed as func- 
D at nominal values of the bulk Reynolds 


number Res for the constant wall temperature boundary condi- 


tions of axial distance x 


tion; the calculations were performed for negative values ol B 
For the constant heat-fl 
values of the initial Nusselt number Nu, 


r/D at nominal values of the initial Reynolds number 


cooling the gas ix boundary condition, 
were computed as fune- 
tions ol 


\ Re 


the gas 


the calculations performed for positive values of B (heating 
were cross-plotted to obtain Nu, for constant @. 
Fig. 8 presents Nu,, with 


computed for air and carbon dioxide 


physical properties varying with temperature 


D, with N pe; and 


as a tunction of 


@ 4s parameters lor the case ol constant heat 


flux (heating). From the afore-mentioned figure it is seen that 


for air the values of Nu, are influenced to a greater extent bv a 
change in @ than are the corresponding values for carbon dioxide, 
It Wiis 


the influence of @ increasing with the Reynolds number 


16 dioxide have almost 


the 


shown in reference that air and carbo 


identical velocit) profiles and boundary layer thicknesses 
calculated values of u,*, t and ts* « 


teria lly different for the Conse qe ntl, 


innot, therefore, be ma- 


two guses , the computed 


value of the Nusselt number is primarily influenced by the varia- 


tion of the thermal conductivity with temperature the variation 
the values 


assumed was (/ = (t/t, of b employed for air 


and carbon dioxide were 0.78 and 1.13, respectively.] 

The effect of b on Nu, for air and carbon dioxide in the case of 
of cooling (8 0 
9. Asin the case of heating, 


heating is also evident in the case as can be 


seen in Fig increasing the intensity 


ol cooling increasing absolute values of 8) causes a much larger 


change in the « omp ited value of Nu, for air than for carbon di- 


oxide 
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Comparison of Theory and Experimental Results 

Fig. 10 compares the local experimental and analytical Nusselt 
numbers for air and carbon dioxide in turbulent flow as a function 
of axial distance and Reynolds number for different values of the 
parameters 8 and @ 
10 


The experimental results presented in Fig 


are those for which the average Reynolds number was in 


closest with 


agreement the nominal values employed in the 
analysis. 

It is seen from Figs. 10(a and b) that for the given 2/D and 
Reynolds number the predicted values* for Nu,° are smaller 
than those obtained experimentally for the case of heating the 
gases, the deviation being larger for carbon dioxide. LOC 
and d@) demonstrate that when the gases were cooled the pre- 
dicted values of Nu, 
mentally, 


anil 


those obtained expe rimentally. 


Figs 


were larger than those obtained experl- 
the difference being somewhat larger at larger values 
of The shape of the predicted curves 18 quite similar to 
It should be noted that the re- 
sults of many investigations are presented in the form h,/h 


Nu/Nu 


ship of the actual values with distance. 


q Ol 
, as a function of axial distance rather than the relation- 
For design purposes, 
however, it is more desirable to have the results presented in 
the form of local values as functions of axial distance, as is done 
here 

Discussion of Disagreement Between Theory and Experiment. ‘The 
the ap- 
those obtaining in the experimental apparatus 


boundary conditions assumed in theory were closely 


proximated b 
The differences between the predicted and experimental results 
are attributed largely to the assumptions needed for computing 
the velocity and temperature profiles (from which Nu° and Np,° 
were calculated 


The following assumptions are those most in question: 


1 The value for @ was assumed to be constant and equal to 
unity 


2 The value employ ed for the diffusivity in the region near 


the wall was assumed to be 


n*uy/» 


n2uy(1 (10 


’ For convenience in discussion, the analytically computed values 
of Nu and Np, are designated with a superscript 
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Fig. 10 Comparison of experimental and predicted local Nusselt numbers for air and carbon dioxide, fully developed initial ve- 
locity profile, uniform initial temperature profile, heating, and cooling 


where the constant 7 0.124 was determined ¢ periment ally fe It has been shown by several Investigators* that the value for a 
idiabatic flow [19 depends upon distance from the wall (however, q@ is relativel 
} The value employed for the diffusivity in the region away constant once the turbulent core has been reached) and upon the 

from the wall was assumed to be that suggested by von Karman — flow conditions (Reynolds number) and Prandtl number. Values 

4 x(du/dy “we 7 or @ that range from approximately 0.2 to greater than 1.5 have 

been reported for Reynolds numbers from 20,000 to 300,000 


where Kk Was assumed equal to 0.38 Consequently, the assumption @ = 1 for all flow conditions 


The constants » and « were determined from experiments with does not agree with indicated experimental results 
adiabatic flow When the flow field is strongly diabatie and the 


fluid properties vary markedly with temperature, it is recog 


Conclusions 


field. Consequently, it is doubtful that adiabatic flow informa- The following conclusions were derived from the experimental 


nized that the velocity field is dependent upon the temperature 


tion can be emploved to calculate heat transfer in strongly dia- and analytical investigations of the local heat-transfer character- 
t ‘ 0 - 0 . 1 y ( \ | f 

batic flow Information is lacking for determining the validity o ‘A discussion of the results of the afore-mentioned investigations 

ind 3% is presented on p. 55, Part II of reference [16 


issumptions 2 
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istics for air and carbon dioxide (Np, 0.7) flowing turbulently 
in smooth round tubes with a fully developed velocity profile 
and a uniform temperature distribution at the cross section where 


energy exchange was initiated. 


1 For the case of heating, values of the one per cent thermal 
entrance length La; ranging from 11 to 27 diameters were deter- 
mined for initial Reynolds numbers ranging from 59,580 to 296,- 
000 for both air and carbon dioxide. For cooling with wall tem- 
perature invariant with axial distance, values of L,,; ranging from 
12 to 26 diameters for air and carbon dioxide were determined for 
bulk Reynolds numbers ranging from 23,000 to 218,700. 

2 The values of L, 


dioxide over the Rey 


for heating and cooling air or carbon 


nolds number range investigated can be ex- 


pressed, with an average deviation of +2 diameters, by means 


the 


of the empirical relationships given by iation (4) with 
constants presented in Table 3. 

3 The predicted values of entrance length L compare favora- 
bly with the value of Li obtained experiment lly 

t The predicted values of local Nusselt numbers for air and 
for carbon dioxide are in moderate agreement with those obtained 
experimentally, heating carbon dioxide 
The predicted values of the local initial Nusselt number Nu,° lor 
heating carbon dioxide 


than the 


except for the case of 


were approximatelh 40 per cent smaller 


corresponding experimental 
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Local and Average Coefficients 


Chemical Engineering Progress, vol. 44, 1948, 


thesis, 


DISCUSSION 


S. Levy° 


The author’s paper is a valuable 


the 
His contribution would be all the more 


iddition to studies of 
thermal! entrance region. 
useful if he sometimes would expand his discussion of results. 


Specific instances include: 


1 Dr. Wolf does not show correlation of the data in terms of 
the standard parameters, Nusselt, Reynolds, and Prandtl num- 
bers. Such correlations would help understand the effects of 
They 


would also make it possible to determine whether the large entry 


fluid property variations for fully established conditions. 
lengths shown in Fig. 3 are not due to changes in fluid viscosity, 
i.e., Reynolds number, in the flow direction. Without such corre- 
lations it is hard to see how the value of h, was obtained since it 
varies along the heated length. Further, comparison of the non- 
dimensional groupings with the data of Humble, et al.,* would 
permit evaluation of the present experimental data and accuracy. 

2 The disagreement between the measured entry lengths for 
heating carbon dioxide and air is difficult to understand (see Fig 
It is 


possible that the difference can be traced to the failure to account 


6). The theory predicts the same value for the two gases. 


Specialist, Heat Transfer and Fluid Dynamics, General Electric 
Schenectady, N. ¥ Assoc. Mem. ASME 
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for variation of the coefficient h, due to the change in local Reyn- 
olds number 

3 The deviation between theory and tests beyond z/D of 10 
to 15 is hard to justify. Deissler’s method is known to give good 
agreement for the case of constant properties and the discrepancy 
needs further amplification. If the analytical method cannot 
account for fluid property variation under fully developed condi- 
tions it must be critically examined before application to the 
thermal entrance region. 

t The variation of the 5 per cent entry length with Reynolds 


number is not mentioned. Is the exponent still the same and 
what happens if the thermal entry lengths are plotted versus exit 
or average Reynolds numbers? 

5 The analysis postulates that the eddy diffusivities for heat 
ind momentum are equal. This is hard to justify in the entrance 
region. Measurements of temperature distribution in the entry 


zone would be very helpful in determining the ratio a 


Author's Closure 


Dr. Levy’s comments raise some very interesting questions 


themselves. It was unfortunate that more discussion could not 
be included, es pec ially regarding the correlation of the aver ge 
heat-transfer results However, since it was necessary to shorten 
the original manus« ript considerably, that would not have been 
possible The following discussion is numbered to correspond 
to Dr. Levy's specific comments 

l rhe correlations obtained for the fully developed heat- 
transfer conditions are presented and discussed in detail in Ref- 
erence [16 Further discussion is presented here under para- 
graph 3. The values for entrance length for heating the gas were 
not determined from Fig. 3, they were determined from the varia- 


T he 


3 are for a constant value of Reynolds 


tion of the « xperiment il heat-transfer coefficient with 2/D 
results presented in Fig 
number, the influence referred to by Dr Levy is due to the varia- 
tior 

In di 


of the 


of the thermal conductivity of the gas with temperature 
termining the entrance length from the experimental values 


heat transfer coefficient (h did achieve a constant value 


the influence of the changing Reynolds number along the length 


of the tube was not accounted for. The results presented in 


Fig. 37 of Reference 16] demonstrated that A for air was a func- 


tion of mass rate of flow per unit area and only slightly influenced 
by wall temperature’ Ah for CO, was a function of mass flow per 
nit area only 


experimental results were corrected, by cross 

data, to a value of Reynolds number equal to the 
ynolds number in the tube it is doubtful that the value 
the 


entrance lengths may be slightly different from those re ported 


iverage 


ol A, would change It may be possible that resulting 
herein 

} With regard to fully developed heat transfer it must be 
d out that Deissler’s method is based on a constant wall 
For the 


the theory pred ts a decrease in Nu at a given 


pointe 


temperature with respect to axial distance case of 


heating (8 > 0 
\ The 


for comparison with his analytical results were obtained in an 


experimental results for air employed by Deissler 


electrically heated tube having very nearly a uniform heat flux 
boundary e« 


dlistunece 


mdition and wall temperature variable with axial 
For the case of heating the gas with the wall tempers 
ture increasing with axial distance (uniform heat flux) the ex 
perime ntal value a of Nu decrease with increasing average wall 
temperature and decrease to a 


increasing 6) for air, smaller 


extent for CO It is unfortunate that Deissler did not have ex- 
tensive experime ntal results, either heating or cooling the gas, 
the 


Fig 


conducted with wall temperature constant with respect lo 


axial distance ll presents the results of the cooling ex- 


periments for air and CO, conducted by the author in the tube 


having an L/D 60 with the wall temperature constant with 
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Fig. 11 Correlation of bulk Nusselt number with bulk Reynolds number 
for air and carbon dioxide; cooling, constant wall temperature 
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Fig. 12 The influence of wall temperature gradient upon the Nusselt 
number for air and carbon dioxide; heating and cooling 


values for fully 
For the ex- 
perimental results of Fig. 11 the values of 8 extended from —0.002 
to —0.10. The influence of 8 on the Nusselt number is small 
compared to that predicted by theory 


pect to axial distance, the results are average 


developed flow hydrodynamically and thermally. 


the theory predicts 
a 40 per cent increase in Nusselt number with a 
in 8 from 0 to —0.1). 


is not immediately apparent. 


tpproximately 
change The reason for this discrepancy 
For the case of cooling the gases 
with positive and negative values of the wall temperature gra- 
dient in the axial direction, see Fig. 12, a negative gradient in- 
creased the Nusselt number, and a positive gradient decreased 
the value of the Nusselt number (at a given Reynolds number). 
The aforementioned results were harmonized by employing an 
extension of the theory due to Rubesin (NACA TN 2345) and 
summarized in a paper by the writer entitled ‘“Heat Transfer 


Transactions of the ASME 





From Nonisothermal Surfaces’’ submitted to ASME, June, 1959, 


but subsequently not accepted for publication. The results 


presented in Figs. 11 and 12 are discussed in detail in Appendix 
D of Reference [16]. 

It is the writer’s conviction that previous failure to correlate 
heating and cooling experiments satisfactorily with a single cor- 
with the 
physical properties evaluated at one reference temperature, may 


relation equation employing dimensionless quantities, 


be attributed for the most part to the different values of wall 
temperature gradients obtained in the experiments analyzed. 


Specifically the cooling data for air and COs could not be corre- 


Journal of Heat Transfer 


lated with the NACA modified Reynolds number method; the 
heating data was correlated with a maximum deviation of + 5 
per cent by employing the term (7',,/7',)° and by evaluating 
the properties at the film temperature, a somewhat better cor- 
relation than that obtained with the modified Reynolds number. 
The experimental data for heating air are in full accord with that 
of NACA Report 1020. 

4 The values of the 5 per cent entrance length are given in 
Table 2 for the range in Reynolds number obtaining in the ex- 


periments. 


The functional relationship of L,; with 


not determined 


Vere Was 
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a siece. ¥ Durbulent Flow in a Circular Tube With 


E. M. SPARROW 


With te tet a Arbitrary Internal Heat Sources and 
Wall Heat Transfer 


An analysis has been carried out to determine the heat-transfer characteristics for turbulent 
flow of a fluid in a circular tube with wall heat transfer. The internal 
heat generation is uniform over the tube cross section, but may vary longitudinally in an 
irbitrary manner. The wall heat transfer may also vary in an arbitrary way in the 
longitudinal direction. The analysis applies along the entire length of the tube, that is 
in thermal entrance as well as fully developed regions. The fluid is assumed to have 


1a 
full 


ner wing 


a 

developed turbulent velocity profile throughout the length of the pipe Numerical 
tre presented for fluids with Prandtl numbers ranging from 0.7 to 100 for 
, embe rs from 50 000 to 500 000 he ¢ vlension of the resul i includ 


e variations is indicated 


Introduction 


results are provided over the teynolds number range from 50,000 


to 500,000 for Prandtl numbers varying from 0.7 to 100 
ips with internal heat generation are Previous analyses [1, 2, 3]! of turbulent pipe flows with internal 


modern technology Such flows may heat generation have been confined to the condition of fully de- 


ieled nuclear reactors and chemical 


' velope d heat transfer, and hence, the results pertain to the portion 
flow passage In many ol these sys- ot pipe beyond the thermal entrance region. Also, these investi 
configuration which is gations were carried out under the restrictions that the internal 
heat gene ration and wall heat transfer are both longitudinal 
lly developed turbulent iniform Laminar pipe flows have, however, lent themselves 


t generation is uniform across more readily to analysis. Results applicable over the entire pips 


manner along the length have been given [4] for the general situation where the 
he transier at ie tube wall is also internal heat generation takes on arbitrary longitudinal ind 


permitted to radial variations, while the wall heat transfer varies arbitrarily 


bitrary longitudinal variation An ex- 


imple of a s tion wh such variations may occur is In a nu- ilong the length. Of course, turbulent flows are of greater im- 
ear reactor he neutron flux varies with positior portance as they are almost always encountered in practice 


The findings of th i lysis can be used to caleulate the heat . 
r characteristics over the entire length of the passage, that General Analytical Considerations 
I 
; rr . . 1! - _ . . on . , 
e therm : well as far downstream. Numerical rhe system to be analyzed is pictured schematically on Fig. 1 
Our attention is directed to the portion ol the tube to the right of 
r ate ) alo ne I 
Pransfer Division and pre — q 0 where the fluid simultaneously experiences an internal heat 
St. Louis, Mo., June 14-18, 1959, pai : 
{ECHANICAL ENGINEERS generation Q(z) and a wall heat transfer g(. The functions Q 
idvanced in papers are to be und q are permitted arbitrary variations with 2 rhe fluid 
3 of their authors and not those of 


ASME Headquarters, March 18, 


, MOV 
ing from left to right, possesses a fullv deve loped turbulent ve- 


‘Numbers in brackets designate References at end 





onstant pres of wall temperature tor case Prandtl! number, v/a 


of wall heat transfer rate of internal heat generation 


series expansion fully developed | temperature 


per unit volume 
ire for case of in distribution for case of uni- 
sources a form wall heat flux area at tube wall 

i heat source fully ke veloped temperature } radial 


local heat-transfer rate per ul if 


co-ordinate measured 


radial pos distribution for case of uni- from tube center line 


form internal heat generation tube radius 


local heat-transfer coefficient dimensionless radial co-ordi- 


in series expansion thermal conductivity nate, r V T./p/¥; Te", al 


mensionless tube radius, 


temperature or cas 


fully developed temperature 
ra heat sources, distribution for case of heat F Vv T)/p/v 
Db @ 


sources dependent on 


teynolds number, 
oefficients in series expansion Nu local Nusselt number, hd /t 


page 


280 NOVEMBER 1959 


Transactions of the ASME 





re) 
or 








— @/Or = q/k at r = ro (specified heat flux 


= 
ae 
—— 





y QO at r = 0 (symmetry) 
r 





at « = 0 (entrance condition 


Fig. 1 Physical model and co-ordinates 
where fg and ¢,.9 are chosen so that 


locity profile which is unchanging with length. The temperature (2a 
icross the entrance section, z = 0, is taken to be uniform at the 
value ft It is desired to determine the heat-transfer character- The two problems as defined by equations (3) and (4) will be 
istics of the fluid at all positions along the tube analyzed separately. 
The starting point of our study is the equation expressing con- We begin by studying the problem of a flow with internal heat 
servation of energy for axially symmetric turbulent flow in a — sources in an insulated tube, that is, / 
pipe. For fully developed hydrodynamic conditions 


be written as 


@ First, consideration is 


, this may given to the situation of uniform heat sources, and a solution is 
obtained which applies both in the entrance and fully developed 
7 > , regions. Then, these results are generalized to the case where 
ri ‘ ; ‘ 
= there are arbitrary variations of the heat sources along the length 
oa Or of the tube. Next, the results for the problem of uniform wall 
heat transfer with no heat sources, originally given in reference 
where a@ and €, represent, respectively, the molecular and eddy ‘ : be: Pith B 
. > a 5|, are reviewed and generalized to apply for arbitrary longi- 
diffusivities for heat lo obtain the energy equation in this form, fies i 
: tudinal variations of heat flux. The linear combination of these 
it has been assumed that the fluid properties are constant, and ' eo ae hed 
: solutions is then made to yield results for the general case of arbi- 
that viscous dissipation and axial heat diffusion are both negligible ; . 

: ; trary longitudinal variations in both internal heat sources and 
compared to heat diffusion in the radial direction . . 
, wall heat flux. Extension of the analysis to the situation of 
Since equation (1) ts linear in the temperature t, we are able to i 

: radially nonuniform heat sources is made in the last section of the 
separate the general problem which includes both internal heat 
report. 
sources and wall heat transfer into two simpler situations in the ot - ; 
. The eddy diffusivity which appears Im equations (56d and ha 
following way a) A problem where there is an internal heat ; C ’ ‘ ; 
ears a word ot IsScUusSs1ION. urrent knowledge ol the relation- 
generation @ in an insulated tube (q U For this case, the os sitels - gid abe 
temperature Is designated as lg by) A problem W here there is a 
wall heat transfer q without internal heat generation Y= UV 


ship between the eddy diffusivity for heat, €,, and the eddy dif- 
fusivity for momentum, €,, is still in a state of uncertainty. How- 

: ait ever, it a ars that for the Prandtl numbers of this study, 0.7 < 
In this instance we denote the temperature | t Then, Irom — si ; “ _ — zs : : 


‘ Pr < 100, the choice of €, = €,, is not unreasonable Previous 
the linearity of the energy equation, the temperature in the com- , ¢ 


his | pr blem i imy ly eneeae tue applic ition of this assumption has led to excellent he it-transfer 
nec yrobiem is s Iply given b 
predictions over this Prandtl number range 


/ f 


oO 


Pe eS | ee ee a oe Uniform Internal Heat Sources in an Insulated Tube 


f may be written as follows: For uniform heating conditions, it is well known that, as the 
flow proceeds farther and farther down the pipe, the heat-transfer 
characteristics approach more and more closely to a limiting con- 
dition This limit is termed the fully developed heat-transfer 
situation and will be denoted by the subseript fd. Taking cog- 
nizance of the fully developed condition, we choose to write the 


temperature distribution tg as the following sum 


,* 





Nomenclature 


absolute temperature; / o a distance from tube entrance ) kinematic \ 
solute wall temperature Ta dimensionless axial co-ordinate, 


absolute bulk fluid tempera x/d 


ture dummy integration variable 


fluid density 
shear stress 


eigentunctior al juation 


l7a 


temperature; &, temperature of radial co-ordinate measuring 
fluid entering tube (a cor distance inward from tube 
stant t*, difference tem wall, 7 } funetion ol Kiven by equation 


perature defined in equatior / dimensionless co-ordinate, 16 


(5); ¢,, bulk fluid tempera , r,/p/v 
| iN f Subscripts: 
ture; ¢,, tube wall tempera molecular diffusivity tor heat, 


ture ] fd denotes fully developed heat-trans- 


pe 


fluid velocity in the x-direction fer condition 


eigenvalues of equation (17a 
a, mean fluid velocity dimensionless total thermal dif- Q denotes situation where Q # 0 and 

= dimensionless velocity in fusivity, (@ + €,)/v q = 0 

direction, u/‘V/ To/p eddy diffusivities for heat and q denotes situation where q # O and 


axial co-ordinate measuring momentum Q 0 
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to” is seen to be a difference temperature which plays an im- 
portant role in the entrance region and approaches zero as z be- 
comes large. We now proceed to solve separately for tg rq and 


, * 


Since energy conservation must 


The Fully Developed Solution. 


everywhere be obeyed, it is necessary that tg ¢4 satisfy equation 


“1 


Ole fd 4 


< 


or pc, 

stant and €, has been replaced by €..0 \ 
fully developed situation for uniform heat 
temperature at all points in the cross section 
mn along the pipe. Expressed in mathematical 


becomes 


Q 


constant = 


I course, satisfy the con- 


Substituting (7b) into (6) and 


inbles, the governing equation for H 


dH 
Pr 


88 counterparts ol u and r as 


nd y represents a dimensionless 


that the 


cessary 
specified. The turbulent velocity 


Ving equations (ret b 


valuated from 


10a 


10b 


26 is taken as the average of equations 
lhe minus-one term appearing on the right 
104) is retained for 26 r > and is de- 


Y, we can proceed with the 


12.8493, evaluated to higher accuracy for the pur- 


tigation, differs slightly from that of reference [6 
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integration of equation (8) in order to find H. An analytical 
solution was not possible and a numerical technique (Kutta- 
tunge method) had to be employed. Solutions were carried out 
on an IBM 653 electronic computer for Reynolds numbers of 
50,000, 100,000, and 500,000 for Prandtl numbers of 0.7, 1, 10, 
and 100. The values of H thus obtained provide, in conjunction 
with equation (7b), the complete solution for the fully developed 
situation. Space considerations preclude tabulation of the solu- 
tions here, but they are available as IBM listings 

A result of practical engineering importance is the wall-to-bulk 
temperature difference for the fully developed situation. First, 
the bulk temperature is introduced by its definition 


‘ Q Qr,? 4(2/d) rr 
I= — 
. pic, k RePr 


Then, substituting into equation (7b) and evaluating H at the tube 
wall(r* =r we find 
(lo, —t P 
20 Se = A(r 12 
Qro?/k 


Using the numerical results for H(1 , the dimensionless wall-to- 
bulk temperature difference has been plotted on Fig. 2 as a func- 
tion of the Reynolds and Prandtl numbers 

We will make further use of the fully develope 1 solution in later 
sections 


The Thermal Entrance Region. the 


temperature distribution in the thermal entrance region is found 


According to equation (5), 


by adding the difference temperature to” to the fully devel ped 


temperature tga. Having determined tg sa, we now turr 
attention to ty* 


To find the governing equation tor ta”, we introduce e ju ition 
(5) into the conservation-of-energy principle (3a Then, taking 
ivccount of the fact that to fd satisfies (6 f 
} 


it follows that tg* must 


where again dimensionless v ariableg have been used It is in- 


teresting to »bserve that the heat source term does not ippear in 
equation 13 

The boundary conditions on tg* may be derived from equations 
3b), in conjunction with (5), giving 


or ‘ Y al rT - 9 l4a 


where we have also us¢€ | equation 75) in finding the last expres- 
s810n. 


We seek a solution for tg” in the form of a product 
= Dor y x 15 


where Y and 2, respectively, depend on z* and alone. By 
substituting into the differential equation (13), we find that 


y = e~ 46x */Re (16 


/ ( do ) ( 282 rt ) ‘ 
d " dr? iar’? 
28? constant 
product solution 


(17a 


the the 


From (14a), the boundary conditions on equa- 


where Re is separation arising from 


tion (17a) are 
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dod/dr+ = 0 at r* =O and (176) 


Equations (17a) and (17b) comprise an eigenvalue problem of 
the Sturm-Liouville type. Solutions are possible 
crete, though infinite, set of 8 values. Hence the 


can be written as 


only for a dis- 
solution for to* 


i > 
Q 3 
= dD. r é ul 
Qro2/k » P| 


n=0 


(15a) 
where (3, are the eigenvalues and @, are th 
and (176 


depend on the 


( orresponding elgen- 
functions of equations (17a 
and 


Since the values of r 
Y in equation (17a Reynolds and Prandtl 
numbers, there will be a different set of eigenvalues associated 
with each choice of these parameters. Utilizing numerical inte- 
gration, the first six eigenvalues and eigenfunctions have been 
found for the same Reynolds and Prandtl numbers previously 
Table 1, while the 
is IBM tabulations 
ire still to be determined 
Applying the boundary condition at the tube entrance as given 
by (145) 


mentioned. The eigenvalues are listed in 
corre sponding eigenfunctions are available 


The coefficients D,, of equation (15a 


we have 


> Dd, = 


n=0 


equation 
LSa 


From this, it follows immediately properties of the 


Sturm-Liouville system that 


It is seen that the computation of the D,, involves an integration 
With the 8 


Liij 


20,000 100,000 


Fig. 2 Fully developed wall-to-bulk temperature difference for uniform 
internal heat source in an insulated tube 


Table 2 Values of the coefficients, E, 


a) Re = 50,000. 


-E 


of the fully developed temperature distribution 


nm? 


,, and D, determined, the difference temperature (g* is knowr 


Irom equation L5a 


, and we may now proceed with the complete 
solution of the problem 

The Complete Temperature Solution. (Combining the results of the 
preceeding sections in accordance with equation (5), we are able 
to write an expression for the temperature distribution which 


applies both in the entrance and fully developed regions, that is, 


1 
).1753x1075 
1242x107 


.1293x1074 


+ 


1328x1075 


4 


-4 


$ z 
RePr d 


f particular practical inte 


Table | 
(a) Re 


‘rest is the 


= 50,000. 


—.- 


2 2 
37 
| 37 

36 
36 


1387 
1380 
1366 
1364 


Re 
65 


68 
67 


ef 

2545 
2538 
2524 
2522 
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67 
) Re 


c 
| 8 


Y 


ize 


LY) 


oOo; Bw 
;}O 


oO 


| 10,670 
10,660 
| 10,660 


oo 
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| 6 ¢ 
4 ‘ 
30 |} 
OS 


Sl 
45 


BS 


12,85 

| 12,790 
12,66 

|} 12,64( 


t 
22 | 
96 


37 


= 500,000. 


T 
| 
| 


r 


| 10,680 | 28,520 | 53,560 
28,490 | 53,490 


28,440 | 53,360 
| 28,430 53,350 


H(r* t bP Do 


difference 


Listing of eigenvalues, 3,,” 


P4 
11,330 
11,220 
10,960 
10,93¢C 


BZ 


20,630 
20,520 
20,270 
20, 240 


6f 


85,810 
85,690 
| 85,440 
| 85,410 


T 


| 


| 


between 


30,190 
30,000 
29,560 
29,510 


66 
125,300 
125,100 
124,700 
124,600 


the 


4 


wall and bulk ten 
uating equation | 


bulk temperature from equation (11 


t 


-6 


19 


Q.u Qt 


Qro?/h 


-E2 


0.5039x10~* 
3616x107 
.3973x1075 
-4263x107° 


(bd) Re 


-Ep 


2454x107" | 0 


-1742107% 
s 
1831x107° 


.1892«1076 


(e) Re 


iperatures at all 


-Es 


| .17S2x107¢ 
.2071x1075 


2377x1076 
1 


100,000. 
-By 
.2228x107¢ 
-8076x10"> 
.8799x1075 


500,000. 


-8, 


+ 
2114x1075 


1489x1075 
Sig 1076 


1537x107 
1537«10 


at the wall (7 r 


Ho 


Ss : 
2, Pid 


n=O 


-9353x10~? 


- 9561416 


stations along the pipe 


vields 


“Es 


+ - - t 
0.2405x10~* | 0.14622x107¢ 
1055x1074 


s 


-1379x1075 
1734x1076 


By -B. 


.4180x107° 
.3056x107° 


-35708x1 


4341x107! 


-E, 


.1158x2075 
8176x1076 


6421x1077 
-8 


kiv 


and substituting t 


-€ 


al- 


he 


1) 


\ convenient rephrasing of this result may be carried out by in- 


troducing the fully developed wall-to-bulk temperature difference 


from (12), giving 


‘Q 


lg u 


where the product 


tos) ta 


D,o,(7 


D> Ew 
, 0 


Hr 


has 


been abbreviated 


18, 
Re d 


| 


a 


by E 


Numerical values of E,, have been listed in Table 2 as a functior 


of the Reynolds and Prandtl numbers, while H(ro* 


directly as the ordinate of Fig 
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may be re 


ad 
2. Using these numerical data, 
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i ~ 24 
a/a@ 


(oe) Pr + 0.7, 1.0 





Prei0 


4 


—— ae a0" 
a/ad 
(b) Pr = 10, 100 


Fig. 3 Wall-to-Lulk temperature variation for uniform internal heat source in an insulated tube 


the o ot te nperature differences given by (20a 
plotted on Figs. 3(a and 6 


has been diameters; while for Re = 50,000 and Pr = 
i function of distance along the 


tube for parametric values of the Reynolds and Prandtl numbers.’ 
The information given on these plots 
Fig. 2, permits rapid evaluation 
lifference 


0.7, the entrance 
length is 21 diameters. A very significant finding is that, es- 
pecially for high Prandtl numbers, these entrance lengths are 
substantially greater than those encountered in turbulent pipe 
flow with uniform wall heat transfer or uniform wall tempers 
ture. These differences in entrance length will be emphasized in 


i later section when the results for uniform wall heat flux are 
ipproach to within 5 per cent of the given 


is 


used in conjunction with 
ol the wall-to-bulk temperature 
it Various stations along the tube 

We can define a thermal entrance length as the heated length 
required for low te to 
fully ce veloped value A dashed line has been drawn on Figs 
Lb 


to facilitate finding the entrance length. So, for ex-  Atbitrary Longitudinal Heat-Source Variations in an 
ample for Re = 50,000 and Pr 100, the entrance length is 15 insulated Tube 
+The curves do not extend all t 


e way to z = 0 because the series We now proceed to generalize these results to apply for arbitrary 
has been truncated 


Ma ane 


variations of the internal heat generation along the tube. For 
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Fig. 4 Representation of an arbitrary longitudinal internal heat-source 
distribution 


the purposes of analysis, suppose that Q(x) can be represented 
as in Fig. 4. Envision for the moment an elementary heating 
process in which there is no internal heat generation up to a 
position %; then, beyond %, there occurs a uniform internal heat- 
ing of magnitude dQ. 
process can be found by applying equation | 


The temperature response to such a 
19). 


(z — &) 
4 
d 


= df 
? RePr 


This elementary heating process can be considered as a differen- 
tial step in the heat-source variation pictured in Fig. 4. By in- 
tegrating the separate effects of these elementary steps, we are 
able to find the temperature distribution at any position cor- 


The in- 


more 


responding to a prescribe upstream variation of Q(z) 
tegral of the right side can be rephr used into a con- 


venient form by integrating by parts. The result of this operation 


Now 


evaiuate 


turning to the wall-to-bulk temperature difference, we 
equation (22) at ind introduce the bulk tem- 


perature by its definition 


finally 


Ziving 


Numerical values of BE, and 8 ire givel 


that the integration can be carried out 


in Tables 1 and 2, so 
either analytically or 
nurrerically) for heat 


any prescribed variation of the internal 


generation Q ry 


Wall Heat Transfer Without internal Heat Generation 
Uniform Wall Heat Flux. 


with uniform wall heat flux has been given in reference [5 


An analysis of turbulent pipe flow 
and 
only the results need be reviewed here 

First, for the fully developed heat-transfer condition, the wall- 
to-bulk temperature difference correspond g to the uniform heat- 
ing. at the wall has been plotted on Fig. 5 for the same Reynolds 
and Prandtl number ranges as heretofore considered. It may be 
noted that the ordinate is proportional to the reciprocal of the 

‘ The results of reference [5] have be 
1 case. 


en extended to cover the Pr = 
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fully developed Nusselt number, which has its usual definition 
head 
Nua = —, ha = — ft _ (25) 
k (tow — beaded 


while the bulk temperature for uniform heat flux is given by 


_ To &(2/d) 


26) 
k RePr ( 


law — tao 
Unpublished experimental air data of Presler and Loeffler [7] also 
appear in Fig. 5 and show excellent agreement with the theory. 
The fluid properties appearing in the experimental correlation 
were evaluated at the film temperature. 
Results for the wall-to-bulk temperature difference which apply 
all along the tube are given in the form of the ratio 


(Ko ~ KA ~ San 


which is plotted on Figs. 6(a and b) as a function of z/d for 
parametric values of the Reynolds and Prandtl numbers. The 
use of Fig. 6, in conjunction with Fig. 5, permits rapid computa- 
tion of those temperature results of greatest practical interest. 
The ordinate of Fig. 6 can also be recognized as the ratio 


Nura/Nu 


where Nu is the local Nusselt number defined by equation (25) 
with subscript fd deleted. 

For this situation of uniform wall heat flux, a thermal entrance 
region can be defined as the length required for (t,,, — t,,) to 
approach within 5 per cent of its fully developed value, or al- 
ternately, as the length required for the Nusselt number to ap- 
proach within 5 per cent of its fully developed value. A dashed 
line has been drawn in Fig. 6 to aid in determining the thermal 
entrance lengths. Comparison of Fig. 6 with Fig. 3, especially 
for the high Prandtl numbers, verifies our previous statement 
about the considerably greater entrance lengths associated with 
the case of internal heat generation. Further, it is interesting to 
observe that the Reynolds number appears to play a much 
greater role in determining the entrance lengths for the heat- 
generation case 


While the important results of reference [5] have been graphed 


| 


EXPERIMENTAL 
DATA (7) 
© T/T, * 1-21-1.24 


Te /Tp + 1.65-1.72 


' -? 
(Taw ad) ig 
Qt, /k 


100,000 
Re 


10,000 


Fig. 5 Fully developed wall-to-bulk temperature difference for uniform 
wall heat flux and no internal heat source 
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Re 
50,000 
—-—-— 100,000 
te 





——-=—— $00,000 


6 24 
a/d 


(b) Pr = 10, 100 


Fig. 6 Wall-to-bulk temperature variation for uniform wall heat flux and no internal heat source 


on Figs. 5 and 6, it will also be useful to give their analytical 
representation as follows 
1f.* xz 


&(2/d Fe Re d 
0 


RePr 


Gi(re* + 


n 


where §,? and F, are listed in Tables 1° and 3, respectively, and 
G(ro*) is read directly from the ordinate of Fig. 5 

Arbitrary Wall Heat Flux. Following a procedure identical to 
that outlined previously, the results for the uniform heat flux 
case as expressed by equation (26) may be generalized to apply to 


* The 8,2 values listed here are slightly refined compared with those 
of ref. [5] 
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arbitrary longitudinal variations g(z). First considering an ele- 


mentary heating process and then integrating the effects of many 
such processes, we find that 


»> _ z - 
ae > FB? [ (Re 
oe Re a" ™ Jy *% 


where 
4 t 
k ReP I — ” 
‘ err 0 


Equation (27) gives the wall temperature at any position for a 


48,7 (x —) 


Re d dz 


Inasmuch as 
the constants 8,? and F,, are available in Tables 1 and 3, the in- 


prescribed upstream variation in wall heat flux. 
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Table 3 Values of the coefficients, F,, 


(a) Re = 50,000, 
-P; Po -Ps = 


0.3507«x1072 
. 2505x1072 
2653x1075 
-2748x10~* 


ae: —.. 

0.1278x107? | 0.9842x1075 
-9686x10"5| .7480x10"5 | .640s«20-3 
3213x1075} 1206x1075 | .110%10-5 

_+1513x10~* | .1614x10°* | .1962«10~¢ } 


0.1848x10~2 
1342x1072 
1846x1075 
-1762x10~* 
(b) Re 
[of [| » | » [| *% | *% | 
7110x1075 0.4400x 1075 

-5136x107> | .s963x1075 | .3277%x107 | 
5887x104 | .4661x10~4 | .4403x10°¢ | 
:6375x10"5| .660910-5 | _.5a36x10~ | 


1 
go 
.1s7sxio"? | .7375x1073 
-1426x1075 | .7971x20-4 
-1455x10~¢ | .8675x10°5 | 


9.5407x1075 


(c) Re = $00,000. 
+ ——y : . 


-Fs 


ai ne 


-Fy | 


“Fs | 


0.1114«1075 | 
-7893x10~4 | 
867% 1075 | 
-876%10°8 | 


-Po 


0.4786x10~5 |0.2590x1075 
3348x1075 | .1618x1075 
-5389*10~* | .1660x10~4 
.340110°5 | .1877x10-5 


1374x1075 
.9715x1074 
1024x1074 
1068x1075 


0.1796«10-3| © 
-1265x1075 
-1312x10~* 
-1333x1075 | 





tegration can be carried out as soon as the variation of q is given. 


Combined Internal Sources and Wall Heat Transfer 


We can now proceed to write the solution for the situation 
where internal heat generation and wall heat transfer occur simul- 
taneously. According to equation (2), results for this combined 
problem are found by adding the contributions due to each of the 
separate problems. 

The most important result of practical interest is the wall tem- 
To 
to the wall and bulk tem- 


perature variation corresponding to given heating conditions. 
compute this, we apply equation (2 
peratures 

eT 


se from the first gives 


fs (30a 


Consideration is first given to the situation where both the in- 
ternal heating and the wall heat flux 


are uniform along the length 


these circumstances, equation (30a) can be recast into a 


form more convenient for rapid calculations as follows: 


tq (306) 


where the subscript fd denotes the fully developed situation 
Numerical values of the first square bracket may be read directly 
from Figs. 3(a and b); while the second square bracket is read 
directly from Figs. 6(a and b). 
differences (lo. w lo. fd ’ found 
and 5, respectively. The bulk temperature for the 


The fully developed tempera- 


ture and (f t»deq are 


from Figs. 2 
situation of combined uniform heat generation and wall heating is 
compute 1 from 


RePr 


So, with this information, equation (305) may be easily evaluated 


to give the wall temperature at various stations along the tube. 


Now, we turn to the situation where the internal heating and 


wall heat transfer vary in an arbitrary way along the length 


Equation (30a) is then applied, and the temperature differences 
(f ‘ 

g to. and t u ty 
equations (24) and (27 


evaluated from 
G2, E,, and 
and 3, 
while the bulk temperature corresponding to the longitudinally 


are, respectively, 


The numerical values of 
F, needed in the computations are listed in Tables 1, 2, 


nonuniform heating is given by 
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z 


tl, = o + q(®)d® + d ((%)dz 


0 0 


(32) 

k RePr 
Once the longitudinal variations are prescribed, equation (30a) 
may then be used to compute the wall temperature at various 
positions. 


Radial Variations of the Internal Heat Generation 

In this paper numerical results have been provided for the 
case of internal heat sources which are uniform across the tube 
cross section. However, the results can be extended to heat 
sources which vary in the radial direction. In this instance, the 
fully developed temperature distribution is given by 


to.ta — too { 


RePr 
Q(r)r dr 


which is of a form similar to equation (7b). The function M(r*) 


must be found by integrating the equation 


2u* rot a 


dM rot Q(r*) 
= 4 + (34) 
RePr ’ dr* Y dr* 


2Pr re* 
= Q(r*)rtdr* 
f, a yrve 


With Q(r*) spécified, this can be solved numerically in the same 
manner and with the same boundary conditions as equation (8). 

The eigenvalues and eigenfunctions arising in the solution of 
the difference temperature (g* remain unchanged, but the D,, in 


Mir?) 


* 
’ 


equation (185) now have to be evaluated using 
of H(r* 


write 


in place 


Hence, for the difference temperature lg we can 


1.* ; 48,47 * 
W ’ ¢ Re 

- = = div d,(r *)e 

2 "ow 

k 0 


f, | M(r*)\r 
D,' = on 
f, rtutd,dr* 


These expressions can be evaluated for any prescribed radial 


where 


"ul o,d) 


heat-source variation, and the extension to longitudinal varia- 
tions proceeds in the same manner as described previously for 
uniform heat sources 
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DISCUSSION 
Stuart W. Churchill® 


At what L/D as one approaches the inlet of the tube do the 
five eigenvalues and coefficients which you tabulate become 
inadequate for computation of the temperature distribution 
and/or heat-transfer rate? 

How do your solutions for heat transfer without internal genera- 
tion compare with those of Sleicher and Tribus?’ 


L. $. Dzung* 


The present paper is a much welcomed supplement to the 
authors’ previous works on less general subjects. The published 
data of eigenvalues and coefficients are certainly very valuable, 
although cases may arise in which higher-order functions than the 
fifth are necessary and it would be more preferable to have a 
method to estimate higher-order terms. 

Such methods are In another connection, 
the writer has solved the problem of prescribed wall heat trans- 


indeed available 


fer without internal heat sources [8],* using Fourier series in- 
It was shown that the 
eigenvalues and the coefficients of eigenfunction expansion can 
be readily obtained analytically if the two radial distribution 
functions of velocity and eddy diffusivity, namely, authors’ ut 


stead of authors’ convolution integral 


and , writer's f = u/i,g = yPr, are known or can be assumed 
Then the eigeniunctions depend only upon two quantities, the 
slope of f at tube wall and an integral 


1 
l= f (2f/g)'/*d(r/ro) 
0 


expressions are 


(37) 


The analytical most suitable for high-order 
eigenfunctions, but low-order terms can also be calculated ac- 
curately with a simple correction 


Kigenvalues for the present problem are, except that of the 


)-« 


where A,,’ is authors’ §,?Pr, / is the integral given in equation 
37), m the order number, and 6’ a correction term negligible for 
high orders 


zeroth order, given by 


(38) 


V \A,'I = r(m 4 


Two qualitative features can be inferred from equation (38) if 


6’ is neglected. Firstly, V4 ,,’ increases linearly with the order 
. / , ° 
m. Secondly, VY A,,'/ is independent upon the type of flow, 
namely, the number Pr or Re. Now the square of the integrand 
of equation (37) can be written as, using authors’ notations, 
u y 


ua v-+ Pr-e,, 


If Re 
expression is negligible compared with Pre,, except near the wall. 
The integral J is therefore roughly inversely proportional to 
the square 
of Pr 


ind Pr are sufficiently high, v in the denominator of this 


root of Pr and hence A,,’/Pr is almost independent 


Both these features are borne out by the authors’ calculations 
From Table 1, the independence of A,,'/Pr = B,? of Pr is at 

* Professor, Chemical Engineering, Department of Chemical and 
Metallurgical Fngineering, University of Michigan, Ann Arbor, Mich 
Mem. ASMI 

(*, A. Sleicher, Jr 


and M. Tribus, “Heat Transfer in a Pipe With 
lurbulent 


Flow and Arbitrary Wall-Temperature 
ASME, vol. 79, 1957, pp. 789-797 

* Brown Boveri & Company, Ltd., Baden, Switzerland. 

* Numbers 8 and 9 in brackets designate additional References at end 
of this discussion 
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once evident. The linear relation is also clear if 8, is plotted 
against the order number. 

This fact seems to indicate that the choice of 8, as eigenvalues 
instead of \,,’ is more fortunate. This may be the case in purely 
turbulent flow, but in a general treatment \,,’ is certainly more 
appropriate since the laminar case is included without modifica- 
tion. 

Authors’ —F,, which is Z,,/A,,’ in writer’s notations, has been 
expressed analytically in a simple equation (71) of reference 
[8] as 


N os 
-~ on Ve) '— Ve 
I Pp m 


where N is a numerical constant, p is the slope of the function f 
at wall, K, is a correction factor approaching unity for high-order 
eigenvalues. 

In case of laminar flow, equation (39) gives results which 
agree perfectly with another publication [9] of the authors. 

Equation (39) suggests that F,, decreases approximately with 
the —*/s power of A,,’. Plotting authors’ values of Z, and F,, 
against 8,* logarithmically, this decrease is indeed observed 
in the case of Z,, although analytic expression of Z, is more 
difficult to work out. The decrease of F,, is less rapid. In fact, 
in some cases F,, even increases with 8,* (Pr = 100). Since the 
series — DF, must converge toward G(ro*+) as explained in refer- 
ence [8], equation (93), this increase is rather doubtful 

If the dependence of J and X,,’ upon Pr as explained in the 
foregoing is taken into account, equation (39) indicates that F,, 
is approximately inversely proportional to the cubie root of Pr 
This is not borne out by the authors’ data. However, since 
other influences have not been considered, this kind of argu- 
ment may not be pushed too far. 


Additional References 


8 L. 8. Dzung, “Heat Transfer in a Round Duct With Sinusoidal 
Heat Flux Distribution,’’ Second International Conference on Atomic 
Energy, Geneva, Switzerland, 1958, Paper P 253. 

9 R. Siegel, E. M. Sparrow, and T. M. Hallman, “Steady Laminar 
Heat Transfer in a Circular Tube With Prescribed Wall Heat Flux.” 
Applied Scientific Research, vol. 7, series A, 1958, pp. 386-392 


S. Levy’? 


The authors have done an excellent job of extending the forced 
convection analyses to the case of internal heat generation. The 
proposed series solutions will undoubtedly yield useful theoretical 
predictions. However, the solutions are only as accurate as the 
assumptions used in their derivations and it may be worth-while 
to dwell upon one of the important approximations made by the 
authors. Siegel and Sparrow assumed that the eddy diffusivity 
for momentum and heat are equal. It is an accepted assumption 
for fully developed flow and thermal conditions, but it remains 
to be justified when the heat generation or temperature vary in 
the flow direction. If von Karman’s similarity theory is used to 
express the two eddy diffusivities, one obtains 


(“*) 
é,=6€, = 
dy 


where / is the similarity scale length and is equal to 


: da / da dt | dt 

dy? / dy dy? / dy 

According to these equations, the temperature 
be proportional to the velocity distribution. 


profile should 
If the two eddy 
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Company, Schenectady, N.Y. Assoc. Mem. ASME 


Transactions of the ASME 





diffusivities are equal such proportionality is not obtained with 
internal heat generation or in the entrance zone. Reference to 
the work of Lin and Shen"! also points out that the foregoing ex- 
pressions for the scale length are based upon fully developed pro- 
files (note the total derivative used in the equation) and are not 
necessarily valid when the temperature and velocity d stributions 
vary in the flow direction. The momentum and heat diffusivities 
are probably not equal at the beginning of the entrance region 
and their difference may even vary with heat generation in the 
flow direction. 


A. B. Metzner‘? 


The authors are to be commended for a well-written and erudite 
analysis of a number of rather difficult heat-transfer problems. 
The following comments, while critical, are not intended to de- 
tract appreciably from the conclusion that the paper is a very 
worth-while contribution to the heat-transfer literature from an 
over-all viewpoint. Specific adverse comments are as follows: 


1 The introduction mentions a number of previous analyses 
that are partially related, but a rather largeamount of high-quality 
work which has been carried out by people interested in viscous 
heat generation in the plastics industry goes unmentioned. In 
particular the publications of Toor,'* Gee and Lyon,'* and Col- 
well,“ appear to be very worth-while contributions whose solu- 
tions are perhaps more complex than the one laminar flow solu- 
tion (reference [4]) singled out for mention in this paper. 

2 The choice of Deissler’s eddy diffusivity variations repre- 
sents what is perhaps the best available choice but nevertheless 
one which is without experimental confirmation in the region near 
the wall. Furthermore, a more recent paper” indicates that 
the eddy diffusivity variations chosen by Deissler predict coeffi- 
cients which are perhaps somewhat low by approximately 15 to 20 
per cent at Prandtl numbers of 100, due to the fact that undue 
weight was placed by Deissler on mass-transfer data which have 
been shown to be subject to some experimental error in the afore- 
mentioned paper. As a result of these uncertainties the results 
predicted in the present analysis must necessarily be used with 
considerable caution until experimental confirmation is obtained. 
The experimental confirmation shown in Fig. 5 is entirely non- 
critical as far as the choice of a set of equations for eddy diffusivity 
is concerned in view of the fact that at a Prandt! number of 0.7 
the variation in the region close to the wall is entirely insignifi- 
cant. 


The present reviewer is puzzled by a statement which appears 
on page 282 of the paper to the effect that ‘‘a result of practical 
engineering importance is the wall-to-bulk temperature differ- 
ence.” The present reviewer does not understand why anyone 
would be interested in this quantity in the case of strictly adia- 
batic flow. Instead it would appear that the entire temperature 
distribution is the interesting result and tabulation of the results 
obtained for the temperature distribution rather than the ab- 
breviated results shown in Fig. 2 might have been the more use 
ful data to show. 


In summary, this reviewer concludes that while there might be 


11C. C. Lin and 8S. F. Shen, “Studies of von Karman’s Similarity 
Theory and Its Extension to Compressible Flows,"’ NACA TN 2541 
and 2542, November, 1951. 

12 University of Delaware, Newark, Del 

13H. L. Toor, ““Energy Equation for Viscous Flow,” [ndustrial and 
Engineering Chemistry, vol. 48, 1956, pp. 922-926; 
ciety of Rheology, vol. 1, 1957, p. 177. 

14 R. E. Gee and J. B. Lyon, ‘“‘Nonisothermal Flow of Viscous Non- 
Newtonian Fluids,” Industrial and Engineering Chemistry, vol. 49, 
1957, pp. 956. 

16 Colwell, Christmas Symposium, 1958 
ciety. 
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also, Trans. So- 


American Chemical So- 


Friend and A. B. Metzner, AIChE Journal, vol. 4 1958, p. 
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some rather specific limitations to the results developed in the 
current paper the authors are nevertheless to be strongly com- 
mended for clearly showing the method of analysis of a difficult 
problem and for developing several specific sets of results in 
sufficient detail that experimental verification can be undertaken 
to determine whether any of these criticisms are serious short- 
comings or of much lesser importance. 


Ralph P. Stein‘’ 


This paper is a very welcome addition to the collection of 
analytical solutions of various convection heat-transfer problems 
presently available in the literature. The combined clarity and 
conciseness of the authors’ presentation is to be commended, with 
one possible exception. This has to do with their explanation of 
how they derive an expression for temperature distribution with 
arbitrary heat-source variation from the uniform heat-source case. 

Apparently, the authors wish to avoid direct use of a formal 
mathematical technique (Duhamel’s theorem), with which Equa- 
tion (22) can be written immediately from Equation (19) but with 
which many of their readers might not be too familiar. Instead, 
they attempt an explanation with more physical significance, an 
approach which is always welcome, if done properly. I am not 
sure that their method of explanation will be easily understood 
by those readers not already familiar with the more formal mathe- 
matical technique. Perhaps clear indication that the integration 
referred to is for the purpose of calculating the temperature dis- 
tribution at any position z, and that the integration is performed 
with respect to the variable 2, which accounts for the heating of 
the fluid as it flows from = 0 to ® = z, would help. 
tally, isn’t dQ or & misplaced in Fig. 4? 


Inciden- 


Authors’ Closure 


The authors extend their sincere thanks to the discussers for 
their interesting remarks and contributions to the paper. 

Turning first to Professor Churchill’s question relating to the 
range of z/d over which the solutions are applicable, we observe 
that in Figs. 3 and 6 the curves have not been drawn all the way 
back tox = 0. 
small value of 2/d as was explained in footnote 3 of the paper. 
This cut-off value represents the smallest x/d for which five 
terms were sufficient to compute the wall temperature results. 
Professor Churchill also makes mention of the analysis of Sleicher 
and Tribus. This work was principally directed to the uni- 
form wall temperature case. The generalization of their results 
to apply to other wall temperature variations cannot be carried 


Rather, all curves have been cut off at some 


out to good numerical accuracy because too few eigenvalues are 
given. Hence, no comparison of their work with the uniform 
wall heat flux case may be made here. 

The remarks submitted by Dr. Dzung represent an important 
adjunct to the paper. It that an immediate 
utility of his findings would be in extrapolating the results of the 


would appear 


present paper to higher eigenvalues. From Dr. Dzung’s dis- 
cussion, it appears that this extrapolation is quite promising for 
the eigenvalues 8,? and the coefficients Z,, but not for the coeffi- 
cients F,. This situation with regard to the F, is not easy to 
understand. We agree that the trend of F, increasing with n 
as exhibited for the higher eigenvalues for high Prandtl numbers 
and low Reynolds number is open to suspicion. But, we re- 
checked the numerical integration by rerunning a few cases with 
a smaller step size, and this did not alter the trend of the coeffi- 
cients. In any event, this behavior is not too important since 
the higher eigenvalues play essentially no role for the high 
Prandtl number results. 


Reynolds number cases, where the step size is not as critical 


For the low Prandtl number-high 


we see no reason why Dr. Dzung’s theory, if correct, should not 


7 Columbia University, New York, N. Y 
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agree with the results of Table 3. It may also be of interest to 
mention an approximate procedure developed by R. D. Cess® 
for finding the higher eigenvalues and eigenfunctions for tur- 
bulent heat Latzko’s 


analysis for the case of flow in a circular tube with uniform wall 


transfer. He applied the method to 


temperature, but presumably other applications are possible. 
rhe remarks of Dr 


need for experimental verification are, 


Levy and Professor Metzner about the 
ol course, quite in or ler 
In this connection, there are some entrance region heat-transfer 
[5] which agree fairly well with predic- 


Addi- 


tional data have recently been taken at the University of Min- 


data shown in reference 


tions based on the current eddy diffusivity expressions. 


nesota'® which also are in fair agreement with the theory. The 


authors feel that Deissler’s eddy diffusivitiy expressions are the 


best that are currently available 


Approximation of the 
Tube Flow 


Heat 


Sciences 


' R. D. Cess 
Transfer in Turbulent 
vol. 26, 1958, pp. 723-724 

9 Redfield W Allen, PhD thesis 
partment of Mechanical Engineering 


Kigenvalues for 
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University of 
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Minnesota, De- 


Now, turning to Professor Metzner’s comments on the wall-to- 
bulk temperature difference, we first note that the flow is not 
strictly adiabatic. Rather, there is a heat source within the 
fluid and hence the bulk temperature varies along the length 
Since the wall temperature is frequently the quantity which one 
wishes to know and since the bulk temperature is usually known, 
the authors felt that 7’, — 7’, was a convenient way of present- 
ing the results. We agree that, for the maximum utility, the 
entire temperature distribution should have been presented. 
However, such a presentation would have required considerably 
more space in the Journal. As noted in the paper more detailed 
information is available as IBM tabulations. 

We thank Dr. Stein for making mention of Duhamel’s theorem. 
This should make the derivation clearer to those who are familiar 
with this mathematical technique. As to Fig. 4, we feel that 
both % and dQ are properly placed (see for example Figure 4.1, 
page 403 of Karman and Biot™). 
von Karman and M. A. Biot, “Mathematical Methods in 
New York, N. Y 


» T 
Engineering,” MeGraw Hill Book Company, In« 
1940. 
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cm. sparrow | Laminar Condensation Heat Transfer 


J. L. GREGG 


. ” 
wa, tows tewere cone, | ON @ Horizontal Cylinder 
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A boundary-layer analysis is made for laminar film condensation on a horizontal 
cylinder. The formulation includes both the inertia forces and energy convection 
terms, which are neglected in Nusselt's simple theory. A similarity transformation, 
valid over most of the cylinder, is found which reduces the partial differential equations 
of the problem (the conservation laws) to ordinary differential equations. Numerical 
solutions of the resulting ordinary differential equations are available for the Prandtl 
number range from 0.003 to 100. Heat-transfer results are presented and discussed. 


Introduction transfer results are presented for the horizontal cylinder over the 
Prandt! number range from 0.003 to 100. 


A THEORY OF LAMINAR FILM CONDENSATION was first 


given by Nusselt in 1916 [1].' Postulating a simple model of the Analysis 


physical phenomenon, he successively analyzed the condensation 


' : The Conservation Equations. The system under consideration 

problem for a variety of geometrical configurations. Attention is eres ae : 
. . , . : is shown in schematic view in Fig. 1. The z-co-ordinate measures 

directed here to the horizontal cylinder. Our aim is to provide a . . . . 

the distanee along the circumference of the evlinder, while y 

modern and more complete analysis of the condensation problem, 7 . 


a , measures the distance normal to the surface. It will be supposed 
bringing it within the framework of boundary-layer theery PI 


The physical model used by Nusselt is as follows: The hydro- that the cylinder is situated in a large body of quiescent, pure vapor 


. which is at its saturation temperature 7's. The surface tem- 
dynamies of the problem were expressed as a simple balance be- ; Bg F = ha 
‘ . , perature ol the cylinder Fis is uniform and Tas < Te A con- 

tween the gravity and the shear forces, inertia effects being dis- : . 
‘ tinuous film of condensate runs downward over the cylinder. 
regarded With regard to the heat transfer, only heat conduc- 
tion was assumed to operate, energy convection not being con- 
sidered. In the boundary-layer formulation, both the inertia 
forces and energy convection are fully accounted for 

The starting point of our study is the boundary layer equations 
appropriate to the horizontal cylinder. It will be shown that 
these partial differential equations can be transformed to ordinary 
differential equations which are valid over a major portion of the 
surface of the cylinder. The transformation can be carried out 


in such a manner that the resulting ordinary differential equations 





coincide with those for condetisation on a vertical flat plate 
Utilizing numerical solutions of the transformed equations, heat- 

'‘ Numbers in brackets designate References at end of paper. 

Contributed by the Heat Transfer Division and presented at the 
Semi-Annual Meeting, St. Louis, Mo., June 14-18, 1959, of THe 
AMERICAN Society OF MECHANICAL ENGINEERS. 

Notre: Statements and opinions advanced in papers are to be 
understood as individual expressions of their authors and not those 
of the Society. Manuscript received at ASME Headquarters, Mareh 
6, 1959. Paper No. 59—SA-13. Fig. | Geometrical configuration and co-ordinates 








Nomenclature 


specific heat at constant pressure a local heat-transfer rate per unit thickness of condensate film 
cylinder diameter area function of X appearing in defini- 
dimensionless velocity function = cylinder radius tion of 7 


mga |, *. ] » rae . . 
acceleration of gravity temperature ae wall te mpera similarity variable, 
latent heat of condensation ture; Tsu, Saturation tempera- (y/r)\ar* y?}' i 


average heat-transfer coefficient ture of vapor AT = Tut — 7 dimensionless temperature, 


(T — Toat)/(7', — Ton 


absolute viscosity 


(Q)/area: AT + = velocity component in x direction 
thermal conductivity = velocity component in y direction 
height of vertical plate dimensionless co-ordinate, 2/7 


. kinematic viscosity 
over-all Nusselt number: Nu r = co-ordinate measuring distance . 


hd/k; Nu, = hL/k along circumference from upper density 

Prandtl number, c,u/h stagnation point shear stress 

static pressure / = co-ordinate measuring radial dis- = function of X appearing in defini- 
over-all heat-transfer rate tance outward from = cylinder tion of F 


unit width surface / = stream function 
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The condensation problem is governed by the basic conserva- 
tion principles: mass, momentum, and energy. The boundary- 
layer form of these laws for constant property, nondissipative, 


1 dod 
Q” A. F , = 
Pr F@ (2 +) 0 


gravity flow over a horizontal cylinder is 


ou Ov 
+ = () 
Or oy 


Ou _— 
+ v = pg sin + pM 
oy r 


O*u 
oy? 


oT ~~) oT ; 
+ F =k (3) 


or oy ov? 


The static pressure gradient Op/Ozr has been neglected in equa- 
tion (2) relative to the gravity force pg sin r/r 

The conservation of mass equation (1) is immediately satisfied 
by defining the stream function W in the usual way, i-e., 


u = ow/dy, v = —ow/dr (4) 


Then, the velocity components u and v appearing in equations (2) 
and (3) are replaced in favor of W, giving 


oy ow oy ow 
= =gsin(r/? rT Vv 
dy Ordy or dy? oy" 
oy oT ow oT k OT 
oy Or Or oy pce, Oy? 


(3a) 


Faced with this formidable pair of partial differential equations, 
we are motivated to seek some means for reducing them to ordi- 
nary differential equations, which are easier to solve 

Reduction to Ordinary Differential Equations. 
experience with boundary laver theory 


Guided by previous 
we propose a new inde- 


pendent variable 7 as follows 
(5a) 


where ¥ is a function of X which is not yet determined. In ad- 
dition, new dependent variables are constructed in the following 


way: 


/ * Fre T- T sat 
Fin)=¥, | | o(X), Hn) = — a 
/ 


(5b) 
v Fa _ T sat 


where @(X) is also undetermined. The variable @ is seen to be a 


dimensionless temperature, while F is related to the velocities of 
the problem by the equations 


v gr* ‘/2 = 
| : | yor 
r v 
v gr* ° do FE d(dy/dX 
r vy? dX y 


The primes denote differentiation with respect to 7 


n” | (5d) 


Inspection 


of equations (5a) through (5d) reveals that the structures of the 


variables 9, F, and @ are those normally found in similar-type 

houndary-laver formulations 
Using these new variables, a formal transformation of equations 
and (3a) 


may be carried out. The end result of this operat ion 


ax * 
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In order that equations (6) and (7) be ordinary differential equa- 
tions, it is required that factors which depend on X be eliminated 
Specifically, it is necessary that 


oy? = asin X 


ee 
o oY y? = bsin X 


ly 
= = ¢ sin X 


oy 

& — (8d) 
dX 
Thus there is seen to be a set 
of four equations for determining the two functions y and @. In 
order that two functions be determined by four equations, it is 
required that the equations are not independent of each other 
but rather, that two of the equations follow from the other two 
Provided that a unique solution can be found, then equations (6 
and (7) will truly be ordinary differential equations and the trans- 
formation defined by equations (5a) and (5b) will have achieved 
its end. 

The problem of finding @ and y from equations (8a) through 
(8d) is precisely the same one faced by Hermann in analyzing free 
convection about a horizontal cylinder [2]. He found that a 
unique solution could be achieved except for small values of X 
However, the boundary layer assumptions are already invalid 
for smal] X for other reasons. 


where a, b, c, and e are constants. 


In view of this, it has been gen- 
erally accepted that the failure to find a unique solution for @ and 
¥ in this region will not introduce errors in the over-all heat trans- 
fer significantly different from those already present in the 
boundary layer analysis. 


The functions ¢ and y as given by Hermann are listed in Table 
12 


Table | and 4 functions 

0° 30° 60° 90° 120° 
0 1.181 2.337 3.430 4.41 
0.760 0.752 0.718 0.664 0.581 


150° 
5.24 
0.458 


165° 180 
5.55 5.84 
0.360 O 


The level of these functions has been adjusted so that the con- 
stants a, b, c, and e of equations (8a) through (8d) are simple 
integers which will lead to a particularly desirable form when in- 
serted in equations (6) and (7). The y function and the product 
Y?, originally presented as figure 17 of reference [2], have been 
reproduced In Fig. 2. 


2 In reference [2], the ¢ function is denoted as f, while the y func- 
tion is denoted as g. 


i ee ae 


| 
2.0} > 





Fig. 2. The functions ) and ¢» (reproduced from Fig. 17, reference [2]) 
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Returning now to equations (6) and (7) and utilizing the @ and 
y functions, there results 


F’"’ + 3FF" — UF’)? +1r=0 (9a) 


0” + 3(Pr)F’ = 0 (9b) 
It is of first importance to observe that these equations coincide 
precisely with the governing equations for film condensation on a 
vertical plate [3]. (This correspondence was assured by adjust- 
ment of the levels of @ and y.) The boundary conditions on F 
and 6, to be discussed, are also the same as for the vertical plate 
problem. As a consequence, the available mathematical solu- 
tions for the F, 6 of the vertical plate are immediately at our dis- 
posal as solutions for the F, @ of the horizontal cylinder. How- 
ever, since different transformations were used, the physical 
variables (such as heat transfer, velocity, and temperature) do not 
foliow the same laws in the two problems. 

The Boundary Conditions. To complete the statement of the 
problem, the boundary conditions must be specified. At the sur- 
face of the cylinder, both velocity components must vanish as a 
consequence of the no-slip condition of viscous flow and of the im- 
permeability of the wall. Also, the condensate immediately ad- 
jacent to the wall must take on the temperature 7. At the 
interface between the liquid and vapor (y = 4), the condition is 
imposed that the vapor offers negligible shearing resistance to the 
condensate flow. This assumption is standard in condensation 
In addition, it will be required that the condensate 
take on the temperature 7',,: at the liquid-vapor interface. The 
formal statement of these conditions is 


theory. 


pou/doy = 0 


T sat 


The thickness of the condensate layer 6, which will vary with 
determined from the 
In terms of the new variables of equations (5a) through 


(5d), these boundary conditions become 


position around the cylinder, must be 
analysis. 


F = 
F’ 5 (10a) 
6 


where ns corresponds to the value of 7 at y 

As a consequence of the transformation (5), our problem has 
been reduced to the task of solving equations (9) subject to the 
boundary conditions (10a). Evidently, these solutions will de- 
pend on the choice of two parameters: Pr and 45. From the 
mathematical standpoint, this is a completely satisfactory situa- 
tion. But, practically speaking, the problem is incomplete, since 
the dimensionless condensate layer thickness ns would not be 
known a priori. So, there remains the need to relate 5 to 
another parameter which contains physical quantities which are 
all known. 

To relate ns (and hence 6) to known physical! quantities, we in- 


voke an over-all energy balance as follows 


é 6 
hty [ pu dy + [ puc,( T sa: _ T dy 
/0 /0 
™ o7 
= h dr (11) 
Oy ‘ 


/ 0 ‘ 
The first term on the left is the energy released as latent heat, 
while the second term is the energy liberated by subcooling of the 


condensate. The right-hand side represents the heat transferred 
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from the condensate to the cylinder over a span of circumference 
from zx = Oto z = zx. In writing equation (11), the standard 
assumption of negligible heat conduction across the liquid-vapor 
interface has been invoked. In terms of the transformed varia- 
bles as given by equation (5), the over-all energy balance becomes 


fF . -tee (11a) 
hig 6'(ns) 

where F(s) and 6’(n5), respectively, denote the values of F and 

d6/dy at n = NB. 

From a solution of equations (9) corresponding to specified 
values of the parameters Pr and 3, the values of F(3) and 6’(m3) 
are available, and the right side of equation (1la) may be evalu- 
ated. So, c,A7'/hi, is determined. For a fixed Prandtl number, 
equations (9) may be solved for a sequence of values of ms, and 
a corresponding set of values of c,A7'/hig may then be computed 
from equation (lla). In other words, for a fixed Prandtl number, 
there is a unique relation between ns and c,AT/hig. So, we may 
regard the solutions as depending on the parameters Pr and 
c,AT/ht,, rather than on Pr and 95. 

Solutions. It has already been remarked that equations (9) and 
the boundary conditions (10a) coincide with those for condensa- 
tion on a vertical flat plate. Further, the relationship between 
c,AT /heg and ng for a fixed Prandtl number, as given by equation 
(11a), is the same as that for the plate. So, immediate applica- 
tion of the mathematical solutions for the plate can be made to 
the horizontal cylinder. 

For the vertical plate, solutions of equations (9) subject to the 
boundary conditions (10a) have been carried out numerically in 
reference [3] for Prandtl numbers of 0.003, 0.008, 0.03, 1.0, 10, 
and 100. The first three of these Prandtl numbers correspond 
to liquid metals, while the last three correspond to ordinary 
liquids. Utilizing these numerical solutions, heat-transfer results 
for the horizontal cylinder will be reported here for values of the 
parameter c,A7'/hyg ranging from essentially 0 to 1 


Results 


Heat Transfer. The quantity of greatest practical interest is 
the heat transferred from the condensate to the cylinder. The 
local heat flux is computed from Fourier’s Law, i.e., 


, ( oT ) 
q = ‘J 
Oy / y= 


In terms of the transformed variables of equation (5), the expres- 


3 ‘ 
- ) ¥' X) (l2a 


where 6’(0) is an abbreviation for (d0/dn)y-o. 


(12) 


sion for gq becomes 
T.) |-0'(0 


The variation of 
q around the circumference is seen to be proportional to the 
variation of the y function. Using Fig. 2, it may be observed that 
q has its maximum value at X = O (upper stagnation point) and 
decreases monotonically with increasing X. 

The over-all heat transfer to the entire cylinder may be found 
by integrating the local values according to the relation 


=2 qdr = 2 ‘ X 
q=2f, q dx rfc aan 


Introducing q from equation (12a) and integrating, there 
tained 


yr  £ 
( - ) p% ie 
- J 


The value of the integral is given by Hermann as 


Q = 2kWTuar — 7, -00)) ( 
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f, y dX = O.6)l0609r 


It is customary to phrase the heat-transfer results in terms of 
un over-all heat-transfer coefficient and an over-all Nusselt aum- 
her defined as follows: 

Q 


= 7 (14 
2ar( Tur — 


With this, the over-all heat-transfer result of equation (13a) be- 


/ hl 6 
Nu, / (™ ) = (.733([—0'(0 (15) 


vy? 


comes 


Since 0'(0) depends upon both Pr and c,A7'/hig, so then does the 


dimensionless heat-transfer representation of equation (15). 
Previous experience with condensation heat-transfer theory sug- 
gests that the dependence on Pr and c,AT he can be sub- 
stantially reduced by rephrasing the dimensionless heat transfer 


in the 


tlorm 


gphial*® ‘ 
Nuz ] PV taf 
kv Tos ] 


Pr 


cypAT /hig |'/* 
= 0.733 [—0'(0)] | pl ] (16 


Utilizing the numerical solutions of equations (9), the right- 
The 
This graph 
contains the complete theoretical prediction of the heat transfer 


hand sides of equations (16) and (lla) can be evaluated 


results thus obtained have been plotted on Fig 5) 


us given by boundary-layer theory 

By inspection of the figure, it is immediately seen that for small 
values of c,AT'/hig, the heat-transfer 
Prandt! number are given by 


gph? V6 
Nu,y = (0.733 t e — 
kw Twi — Ty 


he curves for liquid metals have been extended only to c; AT /hig 
0.1 in consideration of the¥%relatively’small magnitudes of cp/hi 
possessed by these fluids. 


results for fluids of any 


Now, it has been shown in reference [3] that the thickness of the 
condensate film varies monotonically with c, AT /hég 
values of c,A7'/hig correspond to thin films 


So, small 
In very thin con- 
densate films, the effects of both energy convection and inertia 
forces are very small. So, equation (17) describes the heat- 
transfer process in the limiting situation where the energy con- 
vection and inertia play a negligible role. It is interesting to 
observe that Nusselt’s analysis, which deals only with this limit- 
ing situation, leads to a result identical to equation (17) 


that his constant is 1 per cent lower. 


except 


Relative to the prediction of equation (17), the role of energy 
convection would be to increase the heat transfer, while the role 
of the inertia forces would be to decrease the heat transfer. So, 
these two effects are in conflict. With this, we can provide perspec- 
tive for the trends displayed on Fig. 3. For the high Prandtl 
number fluids, where the inertia forces are completely overriden 
by the viseous forces, the effects of energy convection win out and 
the curves rise with increasing film thickness (increasing c,A7’, 
hig For the low Prandtl number fluids, where the heat conduc- 
tivity overrides the energy convection, the effects of the inertia 
forces win out and the curves fall with increasing film thickness. 

Heat-Transfer Results for Small Cylinders. 
condition and a given fluid (i.e., 
theory predicts that 


For a given temperature 


given c,AT7'/hi and Pr), our 


6 = X) = constant 


1 a y( 4 


At any circumferential position X, it is seen that the ratio 6/r will 
increase as the evlinder radius decreases. Since the boundary 
layer equations for the horizontal cylinder are based on the as- 
that 6/r is small, it that 


theory fails for cylinders of very small radius 


sumption follows boundary-layer 

The boundary-laver heat-transfer results may be corrected for 
the effects of small radius in a manner identical to that used by 
Eckert (ref Trans- 
lating his findings into the context of the condensation problem, 


i|, p 165) for the free convection problem 


there is obtained 


Pr-100- 
10a! 
i- 
| 


}— 4 
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Heat-transfer results 
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Nu,’ = ~ 
Inf il +— 
a( r =) 


where Nu,’ is the corrected Nusselt number and Nu, is the pre- 
diction for full-sized cylinders (e.g., Fig. 3 
tions, 2/Nuy< 1, so 


2 2 
Inf l+— ) — Pa 
Nug Nug 


Comparison With Vertical Plate Heat Transfer. 


For normal applica- 


Nug 


For film condensa- 
tion on a vertical plane of height L, the dimensionless heat-trans- 
fer prediction of boundary-layer theory takes the form 


high Vs 
Nu, s 
‘1 LAYCT st — T',) 
( T/h v4 
= (.943[—80'(0 | ra =| (20) 
*r 


where #'(0) is the same function of Pr and cy, A7'/h yy as appears in 
equation (16) and 


AL 
Nu, = J 


So, for a given Pr and cpAT h fo» it follows that 


aphigl.® ‘ 
Nu, / 0.943 = — 
kv T'sat — 7’, 
/ 
Uphigd 
= Nu, / 0.733 — 
/ kv a = 


As a consequence of this correspondence Fig 


o 21) 
ro 


3 can be used 


equally well for the flat plate as for the horizontal cylinder 
It is interesting to inquire about the relationship between plate 
height and cylinder diameter required to produce equal heat- 
transfer coefficients for plate and cylinder. Equating heat-transfer 
coefficients, temperature differences, and fluid properties in equa- 
tion (21), it is found that* 
L = 278d (22 


Using this equivalence, correlations of the heat-transfer coefficient 
for cylinders may also be used for vertical plates and vice versa. 
Variation of Film Thickness. For a particular 
physical situation (i.e., given fluid, A7’ and r), the circumferential 
variation of the film thickness 6 is found from equation (18) to be 


Circumferential 


(23) 


The function ¥\ X) is graphed on Fig. 2 
vides us with the variation of 6. 


and its reciprocal pro- 
> 


Utilizing Fig. 2, 
the condensate layer thickness, starting from 


it is seen that 
a finite value at 
X = 0, grows rather slowly over the upper portion of the cylinder 
and then increases quite rapidly over the lower part of the cylin- 
der. At the lower stagnation point (X = 7 
diction of 


the theoretical pre- 
an infinite film thickness can be interpreted as sig- 
nifying the dropping away of the condensate from the cylinder 
surface 

Comparison With Experimental Heat Transfer. | xperiments on film 
condensation over a horizontal cylinder have been carried out 
primarily for fluids with Prandtl numbers lying above unity. 
Data for these high Prandtl number experiments are summarized 


by MeAdams (ref. [5], p. 340). In general, the conditions of the 


4 According to Nusselt’s simple theory, as quoted by McAdams 


(ref. [5], p. 341) the equality of coefficients is achieved when L = 
2.87 d; a result which is very close to equation (22). 
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tests corresponded to small values of c,AT7'/hig. Since the 
analytical predictions of Nusselt closely correspond to those of 
the boundary-layer formulation for low c, AT’ /hig, the comparison 
of McAdams (table 13-4) between experimental data and Nus- 
selt’s theory also applies here. 

For the low Prandtl number range, experiments on sodium and 
mercury have been carried out by Misra and Bonilla for values of 
cpAT/hig up to 0.03. Their heat-transfer results were only 5 to 
15 per cent of equation (17). While boundary-layer theory pre- 
dicts a lowering of the Nusselt number for these low Prandtl num- 
ber fluids (see Fig. 3), the predicted reduction is much smaller 
than that found by experiment. So, it would appear that the 
effect of the inertia forces is by no means sufficient to explain the 
experimental findings. There still remains the need for further 
experiments to clearly define the departures between the test 
conditions and the analytical model. 
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DISCUSSION 
L. A. Bromley® 


The authors are to be commended for bringing laminar-film 
condensation on a horizontal cylinder ‘‘within the framework of 
boundary-layer theory.’’ I believe it 
compare the result reported here with an approximate one I 


would be of interest to 
derived and corrected after discussions with Rohsenow’ for high 


Prandtl number fluids. 
In the authors’ nomenclature, it may be written: 


| ( 3 =) | 
Nu = (). 7280 Ss hig gph;,d® 
kw To. — T. 


{(. 11 c¢, AT\"* 
10 hig 
the term in braces has the following values as a 
c, AT hte 
c, AT 


p 


function of 


0 01 
ht, 


; } 1.0017 1.017 1.078 1.146 1.260 

From Fig. 3 in the text, it appears that the agreement is good 
for high Prandtl numbers. 

The results for Prandtl number near unity should apply to film 


boiling (within the limits imposed on the derivation) with the 


* Associate Professor of Chemical Engineering, University of 
California, Berkeley, Calif. 
‘L. A. Bromley, Industrial and Engineering Chemistry, vol. 44, 
1952, pp. 2966-2969. 
W. M. Rohsenow, ‘‘Heat Transfer and Temperature Distribu- 
tion in Laminar-Film Condensation,’’ Trans. ASME, vol. 78, 1956, 


pp. 1645-1648. 
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physical properties of liquid and vapor interchanged. Of course, 
the assumption of no drag at the liquid-vapor interface would be 
poor. 

Although I am sure the author is aware of other limitations of 
the derivation as applied to “‘laminar-film” condenser calculations, 
I think it is well to point out some of them to the reader. 


1 Surface tension produces droplets covering appreciable 
areas on the underside of a horizontal tube. 

2 Resistance at the vapor-liquid interface may be present 
either due to ‘accommodation or evaporation coefficient” effects 
or to presence of impurities at the interface. 

3 Ripples may be present. 

t Vapor drag is seldom negligible. 


Authors’ Closure 


The authors extend their thanks to Professor Bromley for his 
interesting discussion. From these comments, it will be possible 


to gain further physical insights into the problem. 
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First of all, it is found that, when the values tabulated by 
Professor Bromley are plotted on Fig. 3, good agreement is 
obtained with the Pr = 100 curve. Now, the Bromley-Rohsenow 
formulation was carried out under the assumption that the 
inertia forces (fluid accelerations) are negligible. This situation 
does occur either for very viscous liquids or for very thin conden- 
sate films. The viscous liquids are characterized by high 
Prandtl numbers; and so it is reasonable that the Bromley 
predictions should agree with the high Prandtl number results of 
the present analysis, where the effects of inertia forces are 
included. Thin films are characterized by small values of 
cyAT'/hy; and according to Fig. 3, the results for all Prandtl 
numbers approach the Nusselt asymptote, which neglects the 
inertia terms. 

The restricting assumptions mentioned at the end of Professor 
Bromley’s discussion are, unfortunately, common to almost all 
condensation analyses. In terms at our present analytical skills 
these problems would be intractable were these simplifications 
not made. 
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The Use of Thin Films for Increasing 
Evaporation and Condensation Rates 
in Process Equipment 


This paper describes an experimental investigation of an evaporating and condensing 
test apparatus in which over-all heat-transfer coefficients as high as 8000 Btu/(hr) 
(sq ft) (deg F) were obtained with water by utilizing thin films both in evaporation 
The films were obtained by wiping on the evaporating surface 
and utilizing surface tension effects on the condensing surface. 


The phenomena 


on both the evaporating and the condensing surfaces are amenable to theory. 


introduction 


= PROBLEM OF DESIGN in heat-transfer equipment 
is generally an economic one influenced by three factors—cost of 
energy consumption, initial cost of equipment, and cost of main- 
tenance. Since both low operating cost and low initial invest- 
ment contribute to an economic optimum, a reduction in over-all 
resistance to heat transfer is a desirable achievement since it al- 
lows either a decrease in the required heat-transfer surface area 
or a reduction in the energy demand by lowering the over-all 
temperature difference across the surface. 

The program undertaken by the authors was aimed at evaluat- 
ing a new concept in heat transfer in which the over-all resistance 
to heat flow is decreased through the combined use of thin films 
The thin films for 
evaporation were achieved by spreading the fluid on the surface 
by means of a slowly rotating wiper 


both for evaporation and condensation. 


Thin films for condensation 
resulted from the utilization of surface tension. The experimental 
program was undertaken during the period from 1955 to 1958. 
\ vapor compressior -tvpe evaporator, using water as a heat- 
transfer fluid, was used to investigate the parameters involved. 
Over-all high as 8000 Btu/(hr) 
sq ft were reached with the experimental apparatus. 


heat-transfer coefficients as 


deg F 


Description of the Heat-Transfer System 


The apparatus used to study the evaporation of thin wiped 
films is shown in Fig. |. It is composed of a heat-transfer surface, 
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Nomenclature 


& vapor compressor, and a rotating wiper. The heat-transfer sur- 
face, upon which water is evaporated and condensed, is a copper 
tube 24 in. long and 3 in. in diameter with a wall thickness of 
0.035 in. The 3-in. copper tube is surrounded by a concentric 
6-in. pyrex glass tube which is wrapped with nichrome ribbon 
and electrically heated to prevent the formation of condensation 
on the inner glass surface. A glass tube was selected to allow 
observation of the condensing surface. 

Inside the center tube is a wiper which rotates about the verti- 
cal axis of the cylinder. Secured to the top of the wiper shaft is a 
water distributor to which raw water is supplied. The distributor 
discharges the water to the inner wall just in advance of the wiper 
so that the liquid flows down the surface in advance of the 
wiper and is spread by the wiper into a thin film over the inner 
surface of the cylinder. As the thin film evaporates, the vapor 
formed is drawn off and transported to the outer annular cham- 
ber by means of a compressor. The compressed vapor, which 
is superheated, is condensed on the outer surface of the center 
tube releasing its heat of vaporization which passes through the 
tube wall and evaporates more water. 


Thin Film Evaporation 


As the wiper rotates inside of. the center tube, it distributes a 
thin film of water for evaporation. Because of the low tempera- 
ture difference across the film, evaporation will occur without 
bubble formation. The agreement of test and theory indicates 
that the heat of vaporization passes through the film by conduc- 
tion alone. The pressure and temperature of the vapor adjacent 
to the evaporating liquid film are slightly lower than the values 
corresponding to saturation at the liquid-vapor interface. For 
evaporation to continue, the vapor leaving the liquid surface 
must be removed from the immediate vicinity. 

The rate of evaporation depends on three factors: 





= thermal conductivity, Btu/(hr 
(sq ft) (deg F per ft 


heat-transfer coefficient, Btu 


hr) (sq ft) (deg F ly ft 
5 = wall thickness, ft 
volume, cu ft 
quantity of heat, Btu 
quantity of matter, lb 


latent heat or 
Btu /Ib 
over-all coefficient of heat trans 
fer, Btu/(hi deg F 
time, hr 


vaporization 


(sq It 
2, ft 
wiper-blade velocity, ft/hr 


saturation temperature, deg | 
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density, |b/cu ft 
length of wiper, ft 
film thickness defined in Fig. 2, 


flow rate, Ib/hi 
area, sq It 


diameter, ft 

speed of rotation, rev/hr 
number of wiper blades 
mass velocity, lb/(hr) (sq ft) 
specific heat, Btu/(Ib) (deg F) 


Subscripts 


linear distances defined in Fig. ) = evaporation 


condensation 
wall 
liquid 
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Fig. 1 
(a) Heat-transfer tube, (b) surface for thin film evaporation, (c) surface for 
thin film condensation, (d) raw water inlet, (e) concentrate outlet, (f) con- 
densate outlet, (g) vapor compressor, (h) vapor compressor drive motor, (i) 
drive motor cooling fan, (j) pyrex glass jacket, (k) raw water distributor, 
(l) wiper, (m) spring-loaded carbon brushes, (n) wiper drive motor 


Cutaway view of experimental heat-transfer apparatus 


1) The initial thickness of the liquid film spread by the wiper 


on the inner wall. This thickness is controlled by the wiper de- 


sign, wiper speed, and wiper load. 
2) The temperature drop across the evaporating liquid film. 


(3) The rate at which the film is renewed on the inner wall 


The simplified analysis of the wiped film which will be pre- 
sented is for a vapor compression-type evaporator similar to that 
shown in Fig. 1. The evaporation system is one in which the 
heat of vaporization is supplied by a vapor condensing on the rear 
side of the evaporating surface. The analysis is based on the 


following assumptions: 
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Fig. 2. Three conditions of wiped thin film evaporation 
1 The heat-transfer surface is sufficiently thin and the wiper 
speed sufficiently low so that effects of heat storage in the tube 
wall can be neglected. 

2 The increase in sensible heat in the evaporating liquid film 
is small compared with the latent heat of vaporization and can 
therefore be neglected. 

3 The heat-transfer coefficient on the condensing surface re- 
mains constant regardless of surface loading and, therefore, for a 
constant over-all temperature difference between the evaporating 
and the condensing surfaces, the heat flow through the evaporat- 
ing liquid film becomes a function only of film thickness 

t The wiper is adjusted in such a manner that the initial 
thickness of the liquid film spread on the inner wall is not a func- 
tion of wiper rpm. 
tube wall is small 
compared with the radial flow and can be neglected. 

6 The 
the boiling point elevation can be neglected 

With these 


ditions of evaporation will occur as shown in Fig. 2. 


5 The circumferential flow of heat in the 


mineral concentration in the water is small so that 


assumptions, it can be envisioned that three con- 
Condition 
a) will take place when the wiper velocity is low so that, at a given 
over-all temperature difference, the wall becomes partially dry 
before the next wiper blade renews the evaporated liquid film. 
Under this condition, only a portion of the available heat-transfer 
surface is utilized. Condition (b) occurs when the wiper velocity 
while keeping the same over-all temperature dif- 
to the point where the entire heat-transfer surface just 
Condition (c) 
the wiper velocity is further increased so that the result is an in- 
crease in the average liquid film thickness between the wiper 
blades 
locity 


is increased, 
ference, 


becomes wet between wiper blades. oceurs when 


The extreme for condition (c) is reached when wiper ve- 
The liquid film thickness between 
wiper blades would then become 


becomes infinite. 
constant and equal to the 
initial thickness of the film spread on the wall. 

Considering first condition (a) in which case the wall is par- 
tially dry between wiper blades, the time required for the water 


film thickness to decrease from its initial thickness J to zero is 
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Fig. 3 Cross section of heat-transfer surface 


equal to the time it takes for the wiper blade to travel a distance 
zr 


For constant wiper blade velocity 
(1) 


where 7, is the time required for a water film to evaporate from a 


thickness & to zero. The evaporation time can be found as 


follow 8s. 
Referring to Fig. 3, the over-all heat-transfer resistance for a 
flat surface is 


The heat transferred through the wall in a time increment dr is 
dQ = UdzLit, — t,)dr (3 


The change in volume of this evaporating water film is 


l dQ 
p hy, 


dV = = Ldxrdd 


Combining equations (2), (3), and (4 


l 6 v t 4 
4 dv = dt 
h. i; phy, 


Integrating equation (5) and inserting the limits 


v = Wo gives 


— ) (Ot — 8 
+ - w,.—- 4 
h, Ky ok ph g 


-quation (6) describes the evaporating water film thickness as a 
function of time, provided the over-all temperature difference and 
the condensing heat-transfer coefficient remain constant. 

The time required for the water film to completely evaporate 


ph,, ( l 6 ) : v2 
T,= —— + VU» + 
t— 6, h. ke 2k; _ 


The wiper velocity can be determined from equations (7 


l 


and 


where z, represents the wetted distance between wiper blades. 
The weight of water spread by one wiper blade is 


w = pLdw (9) 


Combining equations (8) and (9 
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Lz,(t, — t,) 


. 1 i “> 
h 
(5. nt Mea a) 


Equation (10) describes the evaporation rate per wiper blade 
as the wetted distance x, varies between 0 and z. When 2, is 
equal to 2, the evaporation rate reaches a maximum value and 
the film condition corresponds to that shown in Fig. 2(0). 

As the wiper velocity further increases, the average film thick- 
ness increases as shown in Fig. 2(c). For condition (c), the wiper 
velocity can be expressed as 


(11) 


where 7, represents the time required for the film to evaporate 
from J to v;. 


Combining equations (6) and (11) and solving for 7, 
ph 6 vo? — dy? : 
1 = fe ) (% — v,) + : : (12) 
2k 
2k, 


(t, — t,) h, b. 
The time rate of water spread by a single wiper blade is 


w = pL(d> — vi)v 
Combining equations (11), (12), and (13) 


La(t, — t,) 


‘ 1 5 (he + B) 
hy, h T . sa % 
¢ w rg 


Equation (14) gives the rate of evaporation per wiper blade for 
the condition of evaporation shown in Fig. 2(c). 

The average over-all heat-transfer coefficient can be expressed 
as 


wh,,M 


(15) 


The relationship between wiper velocity and rpm for the test 
still is 
v = 7Dn (16) 
In order to evaluate the theory, it was required to have an 
estimate of the initial thickness of the evaporating water film and 
the condensing heat-transfer coefficient. These estimates were 
obtained by simultaneously solving the equations describing the 
wiped film process using test data where wiper rpm was the only 
variable. The calculations showed that the condensing coefficient 
was about 10,000 Btu/(hr) (sq ft) (deg F) and the initial water 
film thickness #» was in the order of 0.0005 in. Assuming that 
these values would remain constant regardless of surface loading, 
equation (6) was used to obtain a set of curves showing water 
film thickness as a function of evaporation time. Equations (6), 
(8), and (16) were then used to find the fraction of the tube wall 
in the test apparatus which remained wet as a function of wiper 
speed, over-all temperature difference, and initial water film 
thickness. Predicted evaporation rates and average over-all 
heat-transfer coefficient as a function of wiper rpm were deter- 
mined from equations (6), (9), (10), (11), (13), (14), (15), and 
(16). 
The curves which result are shown in Figs. 4, 5,6, and 7 Fig. 
4 shows evaporation or distillation rates for the test apparatus 
plotted against wiper rpm for constant temperature difference 
and water film thickness. Note that at constant temperature dif- 
ference the evaporation rate increases linearly with rpm until a 


maximum point is reached. This maximum corresponds to con- 
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dition (b) where the entire tube wall just remains wet. At the 
maximum point the water film resistance to heat flow is a mini- 
mum, the average film thickness being Jo/2. The wall will become 
partially dry if wiper rpm is decreased. As the wiper rpm is in- 
creased, the average water film thickness and the resistance to 
heat flow increase. As the wiper rpm approaches infinity, water 
film thickness approaches a constant value equal to the initial 
thickness of the film spread. 

Fig. 5 shows the evaporation rate plotted against the tempera- 
ture difference for constant wiper rpm and water film thickness. 
The evaporation rate increases nonlinearly with increasing tem- 
perature difference until a maximum is reached. The maximum 
occurs when the entire tube wall just remains wet and the water 
At this 


point, the only way in which the evaporation rate can be in- 


is evaporating from the wall as fast as it is being spread 


creased is to either increase wiper rpm or increase water film thick- 
ness, both combined with an increase in temperature difference. 

Figs. 6 and 7 show the corresponding over-all heat-transfer co- 
efficients plotted against wiper rpm and temperature difference 
Note 
that, in the curves in which wiper rpm is plotted against tempera- 


The area for heat transfer was taken as the wiped area. 


ture difference, the zero temperature difference point corresponds 
to the condition of maximum water film thickness. For example 
7, as the over-all heat-transfer coefficient increases, the 
When the 
maximum coefficient point is reached, the wall is just starting to 
dry 


in Fig 


average water film thickness on the wall is decreasing. 


A further increase in temperature difference has the effect 
of having a greater fraction of the wall dry. As the temperature 
difference approaches infinity, although the distillation rate has 


reached a maximum value, the average over-all heat-transfer co- 
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Fig. 4 Predicted evaporation rates versus wiper rpm for constant 
temperature difference and film thickness 


Curves calculated fort, = 212 deg F, 2 wiper blades each 18 in. long, 3- 


in. diameter copper tube, 5 = 0.035 in., he = 10,000 Btu/(hr) (sq ft) (deg 
F). Refer to Fig. 2 for conditions (a), (b), and (c). 
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Fig.6 Predicted over-all heat-transfer coefficients versus wiper rpm for 
constant temperature difference and film thickness 


Conditions same as in Fig. 4 
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efficient approaches zero at the point where the wall is completel) 
dry. As thinner water film thicknesses are desired which cor- 
respond to higher over-all heat-transfer coefficients, either lower 
temperature differences, higher wiper speeds, or a greater num- 
ber of wiper blades are required to prevent the wall from be- 
coming dry. 

It should be noted that the temperature difference and the 
wiper rpm have no effect on the magnitude of the minimum and 
maximum over-all heat-transfer coefficient. The initial water 
film thickness and the condensing heat-transfer coefficient de- 
termine these values. 


The Influence of Heat Storage in the Heat-Transfer Surface 


The measurements taken on the experimental evaporator were 
found to agree reasonably well with the theory presented, pro- 
vided that: 


(1) The entire wall remains wetted by the evaporating film. 
(2) Part of the wall is dry but the product of wiper rpm times 
the number of wiper blades is low. 


For the condition of a partially dry wall and a high product of 
wiper rpm times number of wiper blades, the effect of heat storage 
in the tube wall was evident. 

Fig. 8 describes the storage phenomenon. The wall is assumed 
to have a resistance considerably lower than the condensing film 
resistance so that temperature differences within the wall can be 
neglected. Also neglected is heat conduction in the axial and the 
circumferential direction and the heat storage capacity of the 
condensing film. 


The lower slanted line shows the temperature gradient after 
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Fig. 5 Predicted evaporation rates versus temperature difference for 
constant wiper rpm and film thickness 


Conditions same as in Fig. 4 
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Fig. 7 Predicted over-all heat-transfer coefficients versus temperature 
difference for constant wiper rpm and film thickness 


Conditions same as in Fig. 4 
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the evaporating film has just dried. The upper slanted line 
shows the gradient after the time interval r. 
The heat transferred to the wall during the time differential dr 


is 
a— Zz 
dQ = ( - ‘) Ah, (t, — t,)dr 
To 


This is equal to the heat stored in the wall 


to — Ly . 
dQ = ( ‘) Adc, p,_dt, 
Lo 


The combination of equations (17) and (18 


(17) 


results in 


h 


c 


- dt 
Pf yO 


For the condition t = 0 at t,, = ty the integration yields 


h 
T 


; 
wo = te — (te — tage whe (20) 


Since the water film resistance is extremely low at the point where 


the film dries, it can be assumed that 


(21) 


The heat stored in the wall is therefore 
o—- Zz rn 
Q = ( ‘) Ac,p,O(t. — t, (, ee ma Mn (22 
Xo 


Lo 
The time 7 and the dry wall fraction ( 


zy — 
- ) can be eliminated 
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Fig. 8 Heat storage in the heat-transfer surface 


Journal of Heat Transfer 


from equation (22) by combining with equations (1), (7), (16) 
and (22). 

The stored heat Q is picked up by the surplus water wave 
which proceeds each wiper and is used for additional evaporation. 
This consideration, which is given in detail in Reference [18],' 
leads to satisfactory agreement between test and theory. 


Condensation and Its Effects on Evaporator Performance 

Tests on the heat-transfer apparatus indicated that extremely 
thin water films could be spread by the wiper on the evaporating 
surface. These films offered such a small resistance to heat flow 
that the maximum over-all heat-transfer coefficient was controlled 
essentially by the maximum condensing coefficient. For this 
reason, a program was undertaken to study the factors influencing 
the rate of condensation. It is known that two types of condensa- 
tion could occur—filmwise and dropwise. 

Filmwise Cond tion on $ 





th Surfaces. The vapor which 
leaves the compressor as shown in Fig. 1 is in a condition of super- 
heat as it surrounds the center tube for condensation. The shape 
of the condensate film which is formed on the smooth surface of 
the vertical tube and the condition of the vapors are shown in Fig. 
9. To determine the average rate of heat transfer under these 
conditions, the analyses of Nusselt [4], Jakob [5], and Rohsenow 
{8] are valid. work includes the effect of 
vapor velocity on condensing film thickness. 

In the analysis presented in this report for thin film evapora- 
tion, the condensing coefficient was assumed constant and equal 
to the average value over the entire length of the tube. Actually 
the average coefficient for condensation on a smooth surface de- 
creases as tube length increases. 


Rohsenow’s recent 


The average theoretical values 
of the condensing coefficient as calculated by the Nusselt method 
are in reasonable close agreement with test. Fig. 10 shows the 
magnitude of the average condensing coefficients that can be ex- 
pected for filmwise condensation on a 3-in. diameter tube 2 ft 
long. 

The initial tests on the compression-type evaporator were 


made with a clean, smooth, copper condensing surface. Since 


1 Numbers in brackets designate References at end of paper. 
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Fig.9 Cross section of the condensate film formed on a smooth surface 
in the experimental evaporator . 
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Fig. 10 Average condensing coefficient versus temperature difference 
for various types of condensation 


Curves (a) dropwise on Tefion, 1 /2-in. diameter tube 3 ft long, ¢. = 212 deg 
F, reference |2|; (b) filmwise on smooth brass, 1/2-in. diameter tube, 3 ft 
long, tc = 212 deg F, reference |19|; (c) filmwise on fluted surface, 2 ft 
= 212 deg F, reference (13); (d) filmwise on fluted surface, 2 ft 
= 86 F, reference (13); (e) calculated from Nusselt theory, 

laminar film condensation on smooth surface, 2 ft long, reference (5) 


the thin films spread for evaporation on the inner tube surface 
offered little resistance to heat flow, the maximum over-all hedt- 
transfer coefficient approached the value of the condensing co- 
efficient. 

Dropwise Condensation on a Smooth Surface. It is well known 
that of the two types of condensation available, drenwise and 
filmwise, the heat transfer under normal operating conditions can 
With this in mind the 
outer surface of the center tube in the test evaporator was highly 


be many times higher with dropwise 


polished and flash chrome plated in order to promote dropwise 
condensation. The test series that followed indicated an increase 
in over-all heat-transfer coefficient with maximum values ap- 
proaching 5000 Btu/(hr) (sq ft) (deg F). This level of dropwise 
coefficients is in general agreement with the values measured by 
other investigators as shown in Fig. 10. 

After about 200 hours of testing on the smooth chrome plated 
tube, the dropwise condensation turned to a mixed condensation 
as shown in Fig. 11, in which part became filmwise. Test data re- 
ported by other investigators indicated that the potential availa- 
ble with dropwise condensation had still not been reached in the 
evaporator setup. Therefore a test program was initiated in 
order to study the effect of dropwise condensation promoters. 
The results obtained with various dropwise condensation pro- 
motors are summarized in Table 1. 

The results of the dropwise condensation studies were dis- 
appointing. The reason is probably the low pressure level at 
which the evaporator was operated in these particular tests. A 
certain amount of uncertainty exists as to the actual effect of the 
Most of the tests re- 
ported in the literature have been made at atmospheric pressure. 
The excellent low pressure tests reported by Gnam [14] are un- 
fortunately clouded by the vapor velocity effects which occurred 
when the vapor pressure was decreased in the test apparatus. 
Kinetic gas theory, however, indicates that with decreasing pres- 
sure level the coefficient of condensation should decrease. 

The conclusion is that our knowledge of the process of drop- 
wise condensation is still rather limited. 


pressure level on dropwise condensation. 


A further extensive 
theoretical and experimental study will be required in order to 
achieve the same understanding as in film condensation. 

Film Condensation on a Fluted Tube. 


302 


If the condensing surface of 


NOVEMBER 1959 


Table 1 Principal results of dropwise promoter tests 


Surface conditions 
Wax on polished 


chrome 


Experimental results 

Appearance of dropwise condensa- 
tion excellent; performance 
extremely low, less than film- 
wise condensation 

Dropwise condensation returned to 
filmwise after five hours 

Dropwise condensation returned 
to mixed condensation after five 
hours 

Dropwise condensation returned 
to mixed condensation after 30 
minutes 

Promoter applied by submerging 
tube in a solution of mercaptan 
and benzene for a period of two 
hours. Excellent dropwise con- 
densation returned to filmwise 
after 30 minutes. 

=xcellent dropwise condensation 
lasted thirty hours 


2) Oleic acid on copper 


Stearic acid on polished 
chrome 


Oleic acid mixed with 
light machine oil on 
polished chrome 

Benzene mercaptan on 
copper 


Mercaptan dissolved in 
1-1-1 trichlorethane 
on copper 

Mercaptan applied di- 
rectly to heated 
copper 


Excellent dropwise condensation; 
returned to filmwise after 30 
minutes 

Promoter applied by submerging 
tube for 12 hours. Excellent 
dropwise condensation for three 
hours. 

The condensate produced was 
cloudy immediately after injec- 
tion. After five minutes of 
operation cloudiness  disap- 
peared. Odor present in the 
condensate. Excellent drop- 
wise condensation; long life 

The condensate produced was 
cloudy immediately after injec- 
tion. After five minutes of 
operation cloudiness disap- 
peared. Odor present in the 
condensate Excellent drop- 
wise condensation; long life 


Dodecanethiol dis- 
solved in benzene on 
copper 


Dodecanethiol dis- 
solved in chlorethan 
on copper 


Mercaptan injected 
into vapor 


a tube is grooved as shown in Fig. 12(a), a condensing heat-trans- 
fer coefficient many times higher than the coefficient for a smooth 
tube can be achieved. The phenomenon was first analyzed and 
reported by Gregorig [13] before the Swiss Physical Society 
in 1952. The effect is due to the pressure gradient produced by 
surface tension in the plane at right angles to the axis of the tube 
Fig. 12(b) shows the profile of the condensing surface as suggested 
by Gregorig. The convexity of the surface at point 1 produces an 
overpressure in the film of the condensate as compared to the 
vapor pressure. The pressure in the condensate is determined b) 
the surface tension and the radius of curvature of the condensate 
film. At point 2, owing to the change in sign of the surface 
curvature, a low pressure is formed in the liquid condensate. 
This produces a pressure gradient in the condensate film in the 
direction from the crest to the base of the groove. The surface of 
the condensate film is therefore sucked into the base of the 
groove. It should be noted that the high rates of condensation 
The base of the 
groove is filled with condensate resulting in a relatively lower heat- 
transfer coefficient in this area. The heat-transfer coefficient 
averaged over the total tube surface is smaller than the crest co- 
efficient but is still considerably greater than that for a smooth 
tube. 

The variation of coefficient with tube length is practically con- 
stant for a given condition of steam and condensate film tempera- 
ture drop. As the temperature difference increases, the con- 
densate in the groove takes on a different profile. 


are confined only to the zone about the crests. 


Each profile 
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Fig. 11 Condensing side of ex- 
perimental apparatus showing 
mixed condensation on an un- 
treated chrome plated surface 
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Fig. 12(a) Fluted tube 
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Fig. 13 Optimized profile dimensions of the fluted surface 


corresponds to a certain condensation rate which is a function of 
tube length. If the temperature difference is increased well be- 
yond the fluted surface design point, a sudden sharp drop in con- 
densing heat-transfer coefficient occurs at the base of the tube as 
the grooves become filled with condensate and the pressure gradi- 
ent is eliminated. When the grooves are not filled, the condensing 


side heat-transfer coefficient is essentially independent of tube 
Fig. 12(b) Profile of the condensing surface used by Cregorig length 
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Based on the concept of using surface tension to maintain a 
A sketch of 


The tube was manu- 


thin film, an improved fluted profile was designed. 
the fluted surface is shown in Fig. 13. 
factured by milling the surface with a profile cutter ground to the 
contour shown 

Evaluation of the fluted surface in the thin film evaporator, 
Fig. 1, produced a maximum over-all heat-transfer coefficient of 
8300 Btu/(hr) (sq ft) (deg F). This was the highest coefficient 
achieved on the test apparatus. The tests indicated that the 
fluted condensing surface had an average condensing coefficient of 
approximately 10,000 Btu/(hr) (sq ft) (deg F) at an over-all tem- 
perature difference of 1.5 deg F. This coefficient was approxi- 
mately four times larger than was achieved with plain filmwise 
condensation on a smooth surface of the same over-all dimensions. 

It is of interest to note that, during test operations with the 
fluted tube, a dropwise promoter was injected and a marked drop 
As the promotor 
washed from the fluted condensing surface with time the high 


in over-all heat-transfer coefficient occurred. 
performance was restored. 


Wiper Design and Its Effect on Evaporator Performance 


The function of the wiper is to spread a thin film of water on the 
This film should be distributed uniformly 
along the entire length of the tube. 


evaporating surlace 
Shearing losses should be 
kept to a minimum so as to maintain high mechanical efficiency. 
The wiper should be constantly supported on a fluid film so as to 
minimize friction and wear. 

During the development program, wiper blades of various de- 
signs were evaluated in the test apparatus. Table 2 presents a 
summary of the designs tested. It was determined early in the 
program that water wettability was the most important wiper de- 
sign criteria 
provide some hydrodynamic or hydrostatic effect, a certain 


Although the design of the wiper is intended to 


amount of contact can exist between the wiper blade and the 
heat-transfer surface which will result in material transfer. If 
the wiper material is not wettable, a hydrophobic film will 
gradually be formed on the heat-transfer surface which in time 
will prevent the formation of a thin wiped film for evaporation 

An analysis of a hydrodynamic slipper bearing type wiper 
similar to the one shown in Fig. 14 indicated that the pressure 
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Fig. 14 Hydrodynamic type wiper in fluted tube 
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generated under the wiper blade depends to a large extent on the 
leading edge boundary conditions. The surplus fluid flowing 
ahead of the blade must also be considered in the analysis. The 
generated hydrodynamic force under the wiper is equal and op- 
posite to a normal force, which can be provided by a spring or by 
centrifugal force or by the combination of both. While the centri- 
fugal force varies as the square of the wiper rpm, for all practical 
Calcula- 
tion of the centrifugal force at low wiper rpm indicated that its 
magnitude was inadequate for the control of film thickness. 
Therefore a spring arrangement as shown in Fig. 15 was used to 
produce the desired pressure. With the spring force large as com- 
pared to the centrifugal force, the thickness of the water film 
spread on the wall is proportional to the square root of the wiper 
Therefore, as wiper rpm increases, the film thickness 
The film thickness can be 


purposes, the spring force can be considered constant. 


velocity. 
will also increase as shown in Fig. 16. 


Fig. 15 Sketch of a spring-loaded carbon wiper 
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Fig. 16 Curve of a hydrodynamic force versus wiper rpm for a hydro- 
dynamic type wiper 
Film thickness curves for water, ff = 120 deg F, B = 0.250 in., carbon 
brush; curve (a)-(a), wipers actuated by centrifugal force; curve (b)-(b), 
wipers actuated by a strong spring force 
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Table 2 Principal results of wiper evaluation tests 
Wiper material 
1/6 in. 


Experimental results 

neoprene sheet Wiper material rubbed on evap- 
orating surface making sur- 
face nonwetting 

1/16 in. neoprene sheet cov- 


Neoprene material took per- 
ered by cotton cloth 


manent set in saturated steam 

at 212 F 

1/,,in. neoprene sheet cov- Short wear life 
ered by cotton cloth and 
spring loaded with 0.005 
inch thick phosphor 
bronze 

0.020 in. thick sheet Kel-F 


Wipers produced thin films pro- 
covered by cotton cloth 


viding excellent heat transfer. 

The wiper life was twenty 

hours. 

0.002 in. brass shim stock Produced extremely thin films. 
High wiper speeds required to 
keep wall wet. 

Wear characteristics good. 
Wiper produced a heavy film. 


0.020 in. thick sheet Kel-F 
covered by two layers of 
cotton cloth and one 
layer coarse weave glass 
cloth. 

0.001 in. thick stainless- 
steel shim stock spring- 
loaded 


Carbon brushes alone 


Produced extremely thin films. 
High wiper speeds required to 
keep wall wet. 

Excellent wear characteristics; 
high wiper speed required to 
produce a thin film 

High wiper speed required to 
produce a thin film 


Carbon brushes __ plus 
weights to increase cen- 
trifugal force 

brushes 
with 


(j) Carbon 
loaded 
spring 


spring- 


Spring force not adequate to 
rubber 


produce a thin film at low 
wiper rpm; rubber took per- 
manent set 

Thin films produced at low 
wiper speeds; long life 


(k) Carbon brushes with 
bronze screening as 
backup spring 


additionally controlled by varying the angle between the wiper 
and the heat-transfer surface. Theoretically, with the centrifugal 
force and the correct spring load in combination, it is possible to 
design a wiper to provide constant film thickness with varying 
rpm. 


Experimental Program 


The experimental evaporator pictured in Fig. | was operated 


for a total of approximately 500 hours. During this time about 50 
combinations of wiper design and condensing surface conditions 
were evaluated. The tests at pressure levels cor- 
responding to saturation temperatures ranging between 120 and 
212 F. 


various wiper-condensing surface combinations 


were made 
Fig. 17 shows several typical sets of test performances for 


Included are 


performance data for the following: 


1 Glass cloth wiper, two blades, smooth chrome plated con- 
densing surface, dropwise condensation, 212 F evaporation. 

2 Glass cloth wiper, two blades, smooth copper tube treated 
for dropwise condensation, 120 F evaporation 

3 Glass cloth wiper, four blades, fluted tube with film con- 
densation, 120 F evaporation 

4 Carbon wiper, four blades, fluted tube with film condensa- 
tion, 120 F evaporation. 


The curves in Fig. 17 show the measured rate of condensation 
per square foot of heat-transfer area plotted against over-all tem- 
perature difference for constant wiperrpm. The over-all tempera- 
ture difference is the difference between the evaporating and 
condensing saturation temperatures. The curves show that the 
condensation rate increases with increasing temperature differ- 
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ence until a maximum over-all heat-transfer condition is reached. 
This point corresponds to the evaporation condition (b) shown in 
Fig. 2. 

Fig. 18 shows a plot of the over-all heat-transfer coefficient 
versus wiper rpm for constant over-all temperature difference. 
At zero wiper rpm some evaporation develops from the liquid 

25 
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Fig. 17 Plot of measured evaporation rate per square foot of wiped area 
versus over-all temperature difference for constant wiper rpm 

Curves (a), (b), and (c) experimental data for cloth wiper, two blades, 
smooth copper condensing surface promoted to give dropwise condensa- 
tion with 1-dodecanethiol and chlorethan, tc = 120 deg F; curve (d) for 
cloth wiper, two blades, stnooth chrome plated condensing surface, tc = 
212 deg F; curve (e) for carbon wiper, four blades, film condensation on 
fluted surface, f- = 120 deg F; curve (f) for cloth wiper, four blades, film 
condensation on fluted surface, ft: = 120 deg F 
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Fig. 18 Measured over-all heat-transfer coefficients versus wiper rpm 

for constant over-all temperature difference 


Experimental data for cloth wiper, two blades, 3-in. diameter copper tube 
chrome plated on condensing surface, 5 = 0.35 in.,t. = 212 deg F, (t. — t) 
= 1.7 deg F. Evaporation occurs at zero wiper rpm as a result of the 
liquid distributed to the wali in front of each wiper blade. With wiper at 
standstill, the liquid spreads as a falling film. 
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which flows down in front of the wiper’s edges. As the wiper is 
Part of the film 
evaporates completely depending on the over-all temperature 
difference and the wiped surface can become partially dry. As 
wiper rpm increases more of the wall remains wet and the average 
over-all heat-transfer coefficient increases until a maximum value 
is reached 


slowly rotated the fluid is spread as a thin film 


\ further increase in wiper speed increases the 
average thickness of the evaporating water film which results in a 
decrease of the distillation rate and the over-all heat-transfer per- 
formance 

Scaling Resistance of Wiped Evaporating Surface. In the search for 
other possible advantages of surface wiping, tests were made with 
natural sea water in order to observe the scaling resistance of 
the evaporator surface. A test was made at 212 F at which sea- 
water evaporators are usually subject to heavy scaling. The test 
No decrease in heat-transfer 
coefficient occurred and a freedom from scale deposits was ob- 
served 


run was extended over 50 hours. 


Conclusions 


The experimental investigation conducted by the authors has 
shown that high over-all heat-transfer coefficients are possible 
through the combined use of thin films for evaporation and 
condensation. The results of this investigation are applicable to 
the design of compact and economically operating process equip- 
ment. Sea-water tests have shown that the slowly rotating wiper 
The 
absence of bubbling in the thin film minimizes brine carry-over 
and as a result a high purity condensate is produced. 


retards the formation of scale on the evaporating surface. 
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DISCUSSION 
Lloyd MacGreagor Trefethen® 


The authors’ findings about the effectiveness of their fluted 
condensation surfaces suggest that such surfaces will be widely 
used in the future. The ability continuously to achieve conden- 
sation coefficients comparable to those of drop condensation is a 
breakthrough of significance to one of the world’s most important 
problems, that of getting enough fresh water. 

I find myself most delighted, however, with the simplicity and 
elegance of the solution, with the “Why didn’t I think of that?’’ 
aspect of those fluted surfaces. Once the principle is pointed 
out, it is clear that fluted surfaces are likely to work very well 
indeed. In fact, one wonders why they were not developed a 
It may well be that the attention given in 
the past to film condensation in two dimensions has impeded our 


decade or two ago. 


thinking; it is obvious, in retrospect, that, by utilizing three 
dimensions, and removing the condensate sidewise, much higher 
coefficients can be obtained, a conclusion that might have been 
suggested years ago by the high coefficients of drop condensation. 
If there is a moral in the authors’ condensation achievements, it 
is probably that huge and unexpected gains can come from de- 
liberate exploitation of all the dimensions that exist in real heat- 
transfer situations as compared to the limited dimensionality of 
the stereotype cases that apparently dominate our thinking. 
The authors do not mention a possibly important aspect of 
the fluted condensing surfaces. The pressure gradient that is 
introduced by the imposed gradient in surface curvature causes 
the free surface as well as the bulk of the film to move sidewise to 
the hollows. Noncondensed gases therefore cannot rest on the 
condensing ridges, but instead will be moved toward the hollows. 
A nonecondensing gas may therefore be less of a poison to conden- 
sation on fluted surfaces than it is to film condensation on flat 
surfaces. 
Their 
statement about drop condensation, ‘“‘Kinetic gas theory indicates 


I should like to question the authors on only one point. 


that with decreasing pressure level the coefficient of condensation 
would seem to presuppose a model for drop 
There appears to be no generally accepted model 


should decrease,’’ 
condensation. 
for drop condensation and indeed the subject seems at the present 
time to be characterized primarily by lack of understanding. 
If there is any validity to a drop-circulation model, which I have 
developed, one would be led to the same conclusion as the authors, 
namely, that the coefficient goes down with decreased vapor 
pressure, but for an entirely different reason; the viscosity of the 
condensed liquid increases. Perhaps the authors could elaborate 
on this point to indicate what properties of the vapor they con- 
sider significant to drop condensation. 


* Professor of Mechanical Engineering, Tufts University, Medford, 
Mass 
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Authors’ Closure 


The authors wish to thank Professor Trefethan for his inter- 
est and comments on this paper. 

Professor Trefethan has pointed out a very interesting aspect of 
the fluted condensing surface: That of its ability to operate in the 
presence of a noncondensable gas. In our continued tests, we 
hope that we will be able to show that this is the case. 

In answer to the question regarding pressure level effect on con- 
densation, we can only state that we agree that there appears to 
be no generally accepted model for drop-wise condensation and 
We do not 


have a definite model which will describe the mechanism of drop- 


that there exists considerable confusion in this field. 


wise condensation. Our statement relating to the kinetic gas 


theory was made with the assumption that the bulk of condensa- 
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We, of 
The kinetic 
We do not 


tion takes place on the surface between the droplets. 
course, have no proof that this is a correct premise. 
theory was used merely to indicate significant trends. 
have available test data to verify this theory. 

The variables which we considered to have influence on drup- 
wise condensation were: The pressure level at which condensation 
occurs, the corresponding saturation temperature, the specific vol- 
ume of liquid and vapor, heat ot vaporization, and the ratio of the 
fraction of the surface which is tree of droplets to the total con- 
densing surface area. 

At the present time we do not intend to delve further into the 
theory of dropwise condensation, since we find that we can achieve 
better results using film condensation on the fluted surface. This 
can be done with less complication than that associated with the 
injection of promoters. 
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Some General Considerations of the 
Heating of Satellites 


Heating of satellites is treated in connection with flight in planetary atmospheres and 
in interplanetary space. Heating during entry into the earth's atmosphere ts generally 
more severe than that encountered during launch or exit, although it can be substanttally 
reduced by employing blunt shapes and small entry angles Decelerations and heating 
rates during entry can be further reduced by employing lifting satellites. These tech- 
niques in combination with known radiation, heat-sink, and ablation cooling techniques 
appear adequate to handle the entry problem at speeds upto escape speed. Return flights 
from distant planets like Saturn and beyond may be characterized by much higher entry 
speeds, and more advanced techniques for coping with the heating problem may be re- 
Heating during entry into the atmospheres of planets neighboring Earth tends 
to be more or less severe than that in the earth's atmosphere depending on whether the 
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quired 


planets are larger or smaller than Earth 

Heating in interplanetary space due to solar radiation becomes severe only tf a satellite 
approaches to within a small fraction of an astronomical unit of the sun, or uf the vehicle 
A cooling problem is encountered in the outer por- 
tion of our solar system, and an auxiliary power source may be required to supplement 
the reduced amount of heat from the sun. The sputtering or erosion of material from 
the surface of a satellite, caused by impact at hypervelocities with small particles in space, 


has an inefficient radiation shield 


may also pose direct or indirect heating problems 


Further research ts required lo assess 


the magnitude of these problems 


Introductory Remarks 


N ORDER TO EXPLORE the heating of satellite ve- 


hicles ntial to understand the environment in which they 


operat f we think of satellites in the broad sense as being 


bodies captured in the sun’s gravitational field, either directly 


by the sun itself, or indirectly by a planet, then this environment 


is our solar system [1}' which is shown schematically in Fig. 1. 


We see that the extent of this system is large indeed, the planet 
Marth being | astronomical unit or almost 100 million miles from 
the sun, while the outermost planet Pluto is about 40 astronomical 
units or 4 billion miles from the sun. Now satellites, as we view 
them, include everything from artificial satellites of Earth and the 
othe 


planets to interplanet iry spacecraft In our general con- 


siderations of heating, therefore, we must contend with an en- 
vironment that includes planetary atmospheres which may be 
relatively dense 


with material, as well as interplanet ry space 


which is relatively void of material but may be relatively dense 
with radiation 

he heating of satellites is largely determined by their motion 
We shall first treat the heating associated with launch or exit of 
shall 


hereafter we shall compare 


4 satellite from the earth’s atmosphere ind then we ittend 


to heating as it arises during entry 
this entry heating proble m with the corre sponding ones likely to 
be encountered during entry into other planetary 


Attention will then be 


itmosphe res 
turned to heating of satellites by solar 
radiation in int rpl inetary space and finally, brief consideration 


will be given to sputtering of surface material from satellites col- 


liding with particles in space. 

! Numbers in brackets designate References 
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or MECHANICAI 

Nott 
understood as individual expressions of their authors and not those of 
the Society. Manuscript received at ASME Headquarters, June 3, 
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MEAN DISTANCE FROM SUN 
ASTRONOMICAL UNITS 


Fig. 1 Solar system 


Heating in Planetary Atmospheres 


Flight of a satellite in a planetary atmosphere tends to be at 


very high velocities because of booster thrust during exit and 
These 


velocities may measure in the tens of thousands of feet per second, 


because of gravitational pull prior to and during entry 


thereby causing intense shearing action between the surface of 
a satellite and the atmosphere. The result is that heat is gen- 
erated at rates and in quantities which may be ample to vaporize 
the entire vehicle unless careful attention is given to its design. 
In order to get a better feeling for this problem and how it can 
be minimized, we have shown atmosphere exit and entry trajec- 


tories schematically in Fig. 2 along with basic heating equations 


[2 Considering first convective heating,? we note that the heat 
rate dH /dt can be reduced by reducing the product of equivalent 
, Which is proportional to Stanton number, 
This 
product is most readily reduced by employing trajectories that 


friction coefficient ( i 


and atmosphere density p at a given flight velocity V. 


*The convective heating equations given are applicable for the 
case of a cool wall, that is, when the wall temperature is small com- 
pared to the recovery temperature. 
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Fig. 2 Trajectories and heating equations 


restrict high flight speeds to high altitudes where there is the 
least amount of air. Thus it is usual to accomplish the exit 
maneuver w-th a locally vertical launch followed by a turn and 
The 


entry maneuver, on the other hand, is characterized by reduced 


acceleration to satellite speed well out in the atmosphere. 


heating rates if it is initiated at a shallow angle to the local hori- 
zontal thereby permitting deceleration at high altitudes, followed 
by a lower speed descent through the denser atmosphere at low 
altitudes. Now heating rates in the stagnation region of the nose 
of a satellite vehicle tend to be high at best because shearing 
action tends to be most severe in th‘'s region. Accordingly, we 
are urged to be especially aware of fundamental factors influencing 
heating rate at the nose, and we see from the second equation that 
the radius of curvature of the surface, o, comes into play through 
Specifically, it is indicated that 
increasing the radius of curvature will markedly reduce stagnation 


its effect on friction coefficient. 


region heating rate dH,/dt, so we anticipate that satellite vehicles 
will tend to have blunt noses. 

We have anticipated the need for reducing heating rates in 
order to alleviate the problem of finding a kind of surface material 
In like 
manner we anticipate a need for reducing total heat transfer in 


or structure, or both, suitable for a satellite vehicle. 


This need 
may become particularly acute in connection with atmosphere 


order to reduce the quantity of material required. 


entry and so we turn attention to the equation for total heating Q 
which holds for this maneuver. It states that total heating is 
proportional to the loss in kinetic energy of the vehicle, and to 
the ratio of friction force to total drag force Cp’S/CpA, where S 
is the surface area over which Cy’ acts and A is the reference area 
for Cp. Since essentially all the kinetic energy is lost in complet- 
ing the entry maneuver, the principal factor of interest to us is 
the ratio of friction force to total drag force. One of the most 
effective ways of reducing this ratio and hence total heating is to 
increase total drag by employing satellites which have highly in- 
clined forward faces in addition to blunt noses. In this manner 
relatively strong shock waves are generated to heat the at- 
mosphere, while relatively weak surface shearing is developed 
to heat the body. 

heat transfer, we will be 
cerned with radiant heat transfer during flight in a planetary 
atmosphere. Radiation to the surface of a satellite from the hot 
gas surrounding it may contribute to heating the vehicle, and 


In addition to convective con- 


radiation from the surface may contribute to cooling the vehicle. 
In either case we are reminded by the radiation equations that 
the rate and amount of such heat transfer are proportional to the 
emissivity and the fourth power of the temperature of the radiat- 
ing medium. To minimize radiation heating rates, therefore, it is 
again attractive to restrict high flight speeds to high altitudes 
where there is less dense atmosphere and hence reduced emission 
of the hot gas surrounding a satellite. Reducing heating rates 
increases the practicality of radiation cooling, of course, because it 
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reduces the surface temperatures required to counterbalance these 
rates by radiation. 

We have, then, a qualitative picture of the heating problem of 
satellites during flight in a planetary atmosphere and some tech- 
niques which may be employed to alleviate this problem. Now 
these techniques did not, in fact, evolve from studies of satellite 
heating, but rather from earlier studies of ballistic missiles [3] 
and other types of hypervelocity vehicles [4] in the earth’s at- 
mosphere. Quantitative studies of the heating of satellites en- 
tering the earth’s atmosphere are of more recent vintage and they 
include those of Kemp and Riddell [5] and Ferri [6], while 
Gazley [7] and Chapman [8] have attended to the heating of 
satellites in other planetary atmospheres as well. The present 
authors have studied these problems [2] with simplified analytic 
techniques suitable for many engineering applications, and it is 
the results of this study which will be reviewed and added to 
herein. Because of the relative importance of heating during 
entry, it will receive the bulk of our attention. 
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Fig. 3. Motion and heating in a launch to circular orbit about earth 

Exit From Earth's Atmosphere. Considering first the exit maneu- 
ver, we have shown in Fig. 3 some calculated motion and heat- 
ing results for an idealized three-stage launch to a near-earth 
satellite orbit. The variations with altitude of velocity and 
flight path angle @ relative to the local horizontal are about as 
we would expect with angle decreasing to zero and velocity in- 
creasing to satellite value at orbital altitude. Launch accelera- 
tions are restricted to 8 g’s or less in this maneuver and orbital 
altitude is somewhat less than 100 miles. Heating rates and hence 
equilibrium surface temperatures for complete radiation cooling 
of the stagnation region of the nose are a maximum at relatively 
These 


rates for a nose of 1 ft radius* are less than 2 Btu per sq ft 


low altitudes and at velocities of about 5000 ft per sec. 


per sec, however, so that surface temperatures do not exceed 
about 2000 deg R which is well within the range of useful strengths 
of molybdenum alloys and other high temperature structural 
materials [9]. In general, then, heating during the exit maneuver 
should be maintainable at levels tolerable to more or less conven- 
tional types of materials and structures. 

Entry Into Earth's Atmosphere. 
during atmosphere entry following an abortive launch and fol- 
lowing normal flight of a satellite through space. 

Entry From Abortive Launch. In the 
launch due to failure of the rocket booster in flight, the vehicle 
may fall back to Earth along an elliptical trajectory reaching 
relatively high velocities and steep flight path angles at entry into 
For launches of the type we have just con- 
18,000 ft per sec tends to yield the 
most critical entry conditions. 


We are concerned with heating 


event of an abortive 


the atmosphere. 
sidered, failure at about 


The conditions in terms of heat- 


3 Bars over o in Fig. 3 and all subsequent symbols or groupings 
thereof shown in this paper indicate that the dimensions of these 
quantities have been omitted; « and CpAa/m have the dimensions 
of feet and cubic feet per slug, respectively. In the discussion of 
these and subsequent heating results, a typical value of unity is as- 
sumed for the barred quantities. 
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4 

65, deg 
Fig. 4 Maximum decelerations and temperatures during entry following 
abortive launch 


ing and decelerations are shown in Fig. 4. 
figure 


On the left of the 


maximum decelerations and temperatures for ballistic 
L/D = 0) are presented as a function of 9,, the 
flight path angle that would have occurred at satellite speed if the 


launch had been successful 


vehicles (i. 


The case 0, = 0 would be a launch 
to a circular satellite orbit, while the larger angles might be used 
The 
maximum deceleration following an abortive launch is about 13 
g's for 6, 


for launch to a parabolic or hyperbolic orbit into space. 


0 and increases quite rapidly to about 20 g's for 6, = 


6 deg. The temperature rise with 6, is small; however, the level 


of temperature is close to 4000 deg Rankine which is approxi- 


mately twice that encountered during a successful launch. Ac- 


cordingly, successful recovery following an abortive launch may 
require the use of ceramics like zirconia for at least the nose of a 
satellite if it is to be radiation cooled, although certainly heat-sink 
and ablation cooling techniques should be considered for this 
region. On the right of the figure we have the variation of 
maximum deceleration and temperature with lift-drag ratio for 
a final launch angle of 0 deg. The temperature decreases only a 
small amount with increasing L/D; however, there is a marked 
reduction in maximum deceleration with increasing L/D (about 
40 per cent for L/D = 1/2). Itis apparent then that the hazards 
of an abortive launch for a manned vehicle may be greatly re- 
duced if a vehicle having a lift-drag ratio of at least 1/2 can be 
utilized. Assuming that we have a vehicle successfully launched, 
we may consider next the motion and heating during atmosphere 
entry from satellite flight. 
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Fig. 5 Effect of entry angle on maximum deceleration 


Entry From Satellite Flight. The maximum deceleration of bal- 
listic vehicles entering at various angles from near-earth satellite 
orbits is shown in Fig. 5. It is indicated that near zero entry angle 
this deceleration is about 8 g’s. Note too that the deceleration is 
only weakly dependent on entry angle for angles up to a degree or 
so; however, for larger angles, the deceleration increases rapidly. 
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Fig. 6 Effect of entry angle on heating 


For manned vehicles it would appear desirable to use the low 
entry angles in order to minimize the decelerations. 

The effect of entry angle on the heating of ballistic vehicles is 
shown in Fig. 6, assuming laminar flow since flight Reynolds 
numbers are generally less than one million [2]. Ratios of heating 
to the heating at zero entrance angle are presented. It is seen 
that increasing the entry angle substantially reduces net heat 
transfer, while the maximum heating rates and hence equilibrium 
temperatures are substantially increased. Thus with a heat-sink 
or ablation-cooled design, higher entry angles may be advanta- 
geous; whereas with a radiation-cooled design it appears de- 
sirable to keep the entry angles as low as possible. At zero entry 
angle, maximum heating rate is about the same as in the abortive 
launch case for 6, = 0, whereas total heating of the stagnation 
region is about twice as great and is equal to 33,000 Btu per sq ft 
for a unit value of the parameter, m/CpAo. Average heating 
rates and average heat transfer per unit area for a complete ve- 
hicle are much less [2, 10]. Thus for example, a vehicle having 
a cone semiapex angle of 45 deg, a nose radius of 1 ft, and a base 
area of 25 sq ft, has an average heating rate and total heat trans- 
fer equal to about '/, the stagnation values for entry at circular 
satellite speed. 
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Fig. 7 Effect of lift-drag ratio on maximum deceleration during entry 





We will consider now the effect of lift on entry from circular 
orbit. The effect of lift-drag ratio on maximum deceleration is 
shown in Fig. 7. It is seen that just as in the abortive launch case 
even very small lift-drag ratios, for example of the order of 1/2, 
markedly reduce maximum deceleration below that experienced 
by a ballistic vehicle. With increasing L/D above 1, decelera- 
tions rapidly approach those of conventional flight. 

The calculated effect of lift-drag ratioon stagnation region heat- 
ing is shown in Fig. 8. Ratios of heat rates and total heat to the 
respective quantities for L/D = 0 are presented as a function of 
lift-drag ratio. We see that increasing lift-drag ratio reduces the 
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Fig. 8 Effect of lift-drag ratio on heating 


earth satellite orbit 
during entry 


heat rates.‘ However, increasing L/D also increases the time of 
flight so that the net effect is to increase the total heat flux. It 
is apparent then that for higher lift-drag ratios, in order to mini- 
mize weight of heat sink or coolant, as much as possible of the 
total heat should be handled by radiation from the satellite. 

We have discussed some factors which influence heating and its 
relation to motion during atmosphere entry at near-earth circular 
satellite speed of about 26,000 ft per sec. Let us consider next 
the problems of atmosphere entry when the entrance velocity is 
escape speed of about 36,000 ft per sec, W hich may be characteris- 
tic of return flights from some of our neighboring planets, stich as 
Venus or Mars, and also our Moon. We will again confine our 
attention to the trajectories characterized by small entrance 
angles and we shall term these grazing trajectories. As in the case 
of entry from near-earth circular orbits, these trajectories lead to 
deceleration in the low density air at relatively high altitudes, 
with the result that radiation heat transfer to the vehicle is 
probably of relatively minor significance [11 Grazing trajec- 
tories can, with accurate guidance, be used to decay a parabolic 
orbit into an elliptical orbit without subjecting a space vehicle to 
excessive decelerations or heating.’ The decay process is shown 
schematically in Fig 9. It is noted that with each graze the 
major axis of the orbit is reduced and is caused to precess 
Several grazes at most are required in order to decay a parabolic 
orbit into a low-altitude elliptical orbit from which the final 
descent can be made with relative ease. Let us examine, in some 
detail, the motion and heating for this type of entry 

The maximum decelerations during a graze for a nonlifting 
vehicle are presented in Fig. 10. In addition, we show exit ve- 
locity as a function of velocity at minimum altitude. For weak 
grazes, the vehicle leaves the atmosphere with little loss in 
velocity and several grazes are required before the final descent 
can be made. For stronger grazes, that is deeper in the at- 
mosphere, the maximum deceleration is larger and if the graze is 
sufficiently strong, the vehicle does not leave the atmosphere but 
merely rises somewhat from the minimum altitude point and then 
falls back in a final descent. For a strong graze, the deceleration 
reaches a peak twice—once in the initial descent, and a maximum 
value in the final descent. It is interesting to note that the 
maximum deceleration during the final descent from a strong 
graze may be less than 8 g’s which is the maximum encountered 
following entry at near-earth satellite speed. Stagnation region 
heating in a graze of a ballistic vehicle is presented in Fig. 11. As 


would be expected, equilibrium suriace temperature and total 


‘ For a given class of configurations, stagnation region heating rates 
tend to be a minimum when m/C Ac is a minimuin and, for constant 
a, this tends to occur in the range of L/D from 1/2 to about 1 [2, 8}. 

Impacting trajectories characterized by larger entry angles may 
be useful for recovering some instrumented packages which can 
withstand very high decelerations provided they are well protected 
perhaps with an ablation heat shield 
vironment [2]. 


from severe entry heating en- 
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Fig. 10 Maximum deceleration during a 
graze 


heat transfer increase as the graze becomes stronger. The 
maximum temperature for a strong graze is about 30 per cent 
higher than for entry at circular satellite speed and the total heat 
transfer is about 50 per cent higher. In order to reduce this tem- 
perature to about the same magnitude as for entry at circular 
satellite speed, the parameter (C,;Aa/m) must be increased by 
a factor of about 8. 

It would appear from these considerations that the use of a 
strong graze for atmosphere entry of ballistic satellites at escape 
speed need not impose prohibitive heating problems nor would the 
decelerations necessarily be greater than 
during entry at circular satellite speed. 


those encountered 


In addition, the time 
duration of the entry for the strong graze is much less than that 
for the several revolutions required in weak grazes. In studying 
the effect of lift on grazing entry, therefore, we will consider 
only the strong graze, and in particular we will consider only the 
use of negative lift. As the vehicle enters the atmosphere the 

® The heat transfer in a strong graze (Vo/V, < about 1.18) includes 
all the heat transferred from entrance to impact. 
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negative lift is used to balance the excess centrifugal force. 
Consequently, the velocity in the graze can be decayed to circular 
satellite speed at much higher altitudes than if positive lift were 
used. This, of course, results in much lower decelerations and 
heating rates than for a strong graze with positive lift. The entry 
procedure is shown in Fig. 12. The vehicle approaches the planet 
on 4a parabolic orbit and enters the atmosphere using negative 
lift. Maximum decelerations and heating occur in the vicinity of 
the minimum altitude point. From here, the vehicle rises and 
approaches a maximum altitude where the velocity has decayed 
to about local circular satellite speed. From this point, positive 
lift is used and the vehicle makes a final descent on a glide trajec- 
tory. 

The calculated effect of lift-drag ratio on maximum decelera- 
tions during a strong graze is presented in Fig. 13. There is a 
marked reduction in maximum decelerations with increasingly 
negative L/D and they are about 50 to 60 per cent lower than 
The effect of lift-drag 
ratio on heating in a strong graze is presented in Fig. 14. The 


those encountered in the final descent 


total stagnation region heating for all lift-drag ratios is approxi- 
mately twice that for entry at circular satellite speed and at 
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Fig. 15 Effect of lift-drag ratio on heating during entry from hyperbolic 
orbit 


L/D 1/2, it is about twice that for L/D = 0. The cor- 
responding equilibrium surface temperatures are about 25 per 
cent higher than for a descent at circular satellite speed, indicat- 
ing that heating rates may be sufficiently high to require the use 
of an ablating material in the stagnation region. Equilibrium 
surface temperatures over the afterportion of a vehicle may be 
sufficiently low that most of the heat may be handled by a radia- 
tion shield with possibly a ceramic coating for the surface 

On a return flight from the vicinity of Saturn and planets even 
further removed from Earth, the entrance velocities may be of the 
In the first 
graze the velocity should at least be reduced to a value slightly 
less than escape speed, with the second graze and final descent 
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being of the types we have just discussed. In this event the maxi- 
mum decelerations in the initial graze would be about four to five 
For lifting ve- 


For example, 


times those for a strong graze at escape speed. 
hicles these decelerations may not be excessive. 
with an L/D of —1/2, the maximum deceleration is about 5 g’s 
and for larger negative L/D it is substantially less. Let us now 
examine stagnation-region heating in entry at superescape speeds, 
We see that both the tempera- 


tures and the total heat are substantially higher than for the 


and this is presented in Fig. 15. 


second graze, temperatures by about 50 per cent and total heat 
by about 75 per cent. The problems of entry at superescape 
speeds are indeed severe, and it is probable that more advanced 
cooling techniques will be required to cope with these problems 

Entry Into Other Planetary Atmospheres. [n closing our discussion 
of atmosphere entry it is appropriate to consider satellites enter- 
ing the atmospheres of planets other than Earth. 
relative decelerations and heating encountered during entry into 
These 
|, and the more recent and 
indicated that the 
problem of landing on Venus is, from the standpoint of decelera- 


Estimates of 


atmospheres of various planets are summarized in Fig. 16. 
results agree with those of Gazley [7 
extensive studies of Chapman [8]. It is 
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Fig. 16 Relative deceleration and heating 


Fig. 14 Effect of lift-drag ratio on heating during 


tions, comparable to that for landing on Earth, while that of 
landing on Mars is, in this respect, substantially less severe. 
Entry into Jupiter’s atmosphere, however, is much more severe 
and perhaps should not be attempted with other than instru- 
mented vehicles. More or less the same situation is predicted 
for heating. However, precise estimates of the heating associated 
with entering the atmosphere of any planet other than Earth 
must await the accurate determination of the properties of these 
It is useful to note, nevertheless, that the problem 
of entering the earth’s atmosphere appears, from the motion and 


atmospheres. 


heating point of view, to be equally if not more severe than that of 
entering the atmosphere of Venus or Mars. Accordingly, designs 
suitable for entering our atmosphere may also be adequate for 
entering these other atmospheres. A qualification of this state- 
ment should be made in the case of entering the atmosphere of 
In this event 
radiation heat transfer to a satellite may be far more severe than 


Venus which may be largely carbon dioxide [12]. 


in the earth's atmosphere, depending on the emissivity of COs, 
under entry conditions. 


Heating in Interplanetary Space 

Heating by Solar Radiation. While radiant heating is frequently, 
if not commonly of secondary importance during atmosphere en- 
try, it is of primary importance during space flight. The principal 
For 
most of our purposes the sun may be considered a point source of 
radiant energy with the heating rate inversely proportional to the 
square of the distance from the sun. 


source of radiant energy in this case is, of course, the sun. 


The heating of unit area of 
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Fig. 17 Radiant heating by sun 


i satellite which is always considered to be facing the sun has been 
calculated’? for one revolution around the sun The results of 
these calculations for vehicles launched from Earth are shown in 
Fig. 17 which includes radiation equilibrium temperature, total 
heat transferred, and the radiation heat rate as a function of peri- 


helion distance from the sun, The total heat increases quite 
rapidly as the perihelion distance is decreased and the amount is 
generally so large that structural designs based on the heat-sink 
principle with no radiation cooling do not ippeal feasible. For 


approaches as close as astronomical unit the heat rates and 
hence the equilibrium surlace temperatures do not exceed those 
encountered by vehicles entering the earth's atmosphere at circu- 
lur satellite speed It is indicated, therefore, that with radiation 


shields of ceramic design a satellite could eathe 


exceedingly 


The heat 


follows 


ng equations are as 


Radiation rate 0.119 , Btu/tt 
remperature, Type’ + 706/ 4 les 
lotal heat, H/re\ S88 & 106/(h 


where 


ION BEAM APPARATUS 


close approaches to the sun and, if properly instrumented, such a 
sun probe might provide highly valuable data. For satellite flights 
to regions far more removed from the sun than Earth, our “‘so- 
called’? heating problem will in fact become a cooling problem. 
For example, in the vinicity of Jupiter, the solar heating rate 
is only 0.0042 Btu/ft?-see which is approximately 3 per cent of the 
value at the earth and the radiation equilibrium temperature 
(Tree /*) is —160 deg F. 
sun may have to be supplemented with energy from an auxiliary 


The reduced amount of heat from the 


power source. 

Sputtering by Hypervelocity Particles. Another source of heating of 
satellites in space may be due to hypervelocity collisions with 
particles like micrometeorites or gas molecules. Surface erosion 
due to hypervelocity impact with gas molecules has been termed 
“sputtering.’’ In general, the rate at which material is sputtered 
or eroded from the surface and the heating due to impact with 
these particles are small. However, sputtering may limit the use- 
ful lifetime of satellites, and could conceivably roughen the sur- 
face so as to induce turbulent flow during the critical heating 
period of atmosphere entry, thereby making entry heating even 
more severe 

Experiments to determine the nature of sputtering of various 
materials are being conducted at the Ames Research Center. 
The equipment used is the Ion Beam Apparatus shown in Fig. 18. 
The apparatus consists of an ion source, accelerator, analyzer, and 
focusing system. The source yields molecular and atomic oxygen 
and nitrogen ions with intensities in the milliampere regime. The 
accelerating system provides a continuously variable range of 
a few electron-volts to 9 k electron-volts, corre- 
sponding to speeds of about 30,000 to 800,000 ft per sec for ni- 
trogen molecules. 


energies from 


The flat faced aluminum model shown at the 
right was about 1 centimeter in diameter and was exposed to an 
ionized nitrogen stream for 3 hr at an energy level of 3.0 kev cor- 
responding to a stream velocity of about 460,000 ft per sec. A 
total charge of about About 


grams of material were removed and evidence of surface roughness 


1 coulomb was used. 150 micro- 


is shown in the view on the lower right. The velocity of the ex- 
periment corresponds to that of a vehicle at a perihelion distance of 
about 0.09 astronomical unit or about 8 million miles from the sun. 
Calculations indicate that 


sputtering, even at these velocities, 


would be smail provided the surface temperature 18 moderately 


BEFORE 


Fig. 18 Sputtering at hypervelocities 
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low. Surface temperatures at this distance from the sun may, 
2500 deg R. At elevated temperatures 


materials lose their strength and probably their resistance to 


however, reach about 


sputtering, so the effects of sputtering may become more severe 
near the sun 


Concluding Remarks 


Recapitulating our general considerations of satellite heating, 
we observe that they can logically be divided into two parts, one 
dealing with flight in planetary atmospheres and the other dealing 
with flight in interplanetary space. Convective heating during 
launch or exit from the earth’s atmosphere appears to be rela- 
tively modest provided high flight speeds are restricted to high 
altitudes where air density is low. Reduced entry angles sub- 
stantially reduce heating rates during atmosphere entry, although 
the rates remain severe at best by comparison to those en- 
countered during exit. Lifting satellites experience much lower 
heating rates and decelerations, which may make them especially 
Heating during atmosphere 
entry from abortive launches, near-earth circular orbits, and near- 


attractive for manned applications. 


planetary orbits appears tractable if blunt shapes and known 
radiation, heat-sink, or ablation cooling techniques are used. How- 
ever, more advanced cooling techniques may be required to cope 
with the heating during entry following a space trip from the 
outer planets. Entry into the atmosphere of Venus appears to 
be comparable, from the heating and deceleration point of view, 
to entry into the earth’s atmosphere, while entry into the at- 
mosphere of a smaller planet like Mars appears to be much less 
severe, and entry into the atmosphere of a larger planet like 
Jupiter appears to be much more severe, 

Heating by solar radiation of a satellite in interplanetary space 
becomes severe only if the vehicle approaches to within a small 
fraction of an astronomical unit of the sun, or if it is equipped 
Interplanetary flight to the 
outer portions of our solar system will pose a cooling rather than a 


with an inefficient radiation shield. 
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heating problem, and an auxiliary power source will probably be 
required to supplement the reduced amount of heat received from 
the sun. Sputtering or erosion of material from the surface of a 
satellite impacting at hypervelocities with small particles in space 
may also pose problems in the form of reduced lifetimes of the 
vehicles or increased entry heating because of surface roughness 
effects. 
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The Penetration of Planetary Atmospheres 
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The problems of penetrating a planetary atmosphere are reviewed and techniques for 
penetration are discussed. 
Gradual entry into the atmosphere with relatively low deceleration loads and heating 
rates low enough so that the heat may be rejected by thermal radiation from the surface; 


These techniques are resolved into two general classes: (1) 


and (2) direct entry with higher deceleration loads and higher heating rates with the heat 
being absorbed by the body surface. 


Introduction 


HE PROBLEM OF atmospheric penetration arises for 
any vehicle which approaches a planetary atmosphere and for 
which physical recovery or survival is desired at the planetary 
surface. This includes a gamut of possibilities from simple 
sounding rockets to manned vehicles returning from inter- 
planetary trips. Several types of atmospheric entry are il- 
lustrated in Fig. 1. 
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Fig. 1 Various types of atmospheric entry 

The approach to a planet from space is presumed to involve the 
classical technique of first getting the vehicle into the planet’s 
solar orbit and then “falling’’ into the planet. The approach 
velocity is, then, essentially the escape velocity of thé planet 
corresponding to a parabolic orbit. (In the vicinity of the earth, 
the escape velocity is about 37,000 ft/sec; a highly eccentric el- 
liptical orbit, such as a circumlunar orbit, involves a velocity 
almost as high—about 35,000 ft/sec.) It is expected that rather 
precise control of velocity and direction would be required to ob- 
tain a hit or a near-miss. A direct hit would involve an entry 
path similar to the ballistic rocket—but with a higher initial ve- 
locity. A more gradual descent can be accomplished by either a 
tangential shallow entry or by first maneuvering into a satellite 
orbit and then descending. 

Although the heating and deceleration accompanying at- 
mospheric entry bring about severe design problems, the presence 
of a planetary atmosphere is advantageous in that it acts as a 
cushion to reduce a space vehicle’s velocity to a safe landing 
speed. Otherwise, one would be forced to the costly expedient of 
rocket braking. Landing on the moon, which has effectively no 
atmosphere, requires an appreciable portion of the vehicle’s 
weight to be allocated to a braking rocket. Landing on a planet 
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with an atmosphere requires a clever approach shot, a surface- 
cooling scheme, and/or a clever aerodynamic design. 


Characteristics of Atmospheric Penetration 


A body approaching a planetary atmosphere possesses a large 
amount of energy—potential energy because of its position above 
the planet’s surface and kinetic energy because of its velocity. 
In most cases, the kinetic form predominates. A satellite in a 
circular orbit at 200 mi altitude possesses a kinetic energy of about 
13,000 Btu/lb. An interplanetary vehicle approaching the earth 
possesses a kinetic energy of about 26,000 Btu/lb. That these en- 
ergies are not negligible is illustrated by comparison with the 
energy required to vaporize water, 1000 Btu/lb, and to vaporize 
carbon, 25,000 Btu/lb. The latter is among the highest heats 
of vaporization of known materials. 

One of the problems of atmospheric penetration is to convert 
this energy into a form that will not prove disastrous to the ve- 
hicle, either during penetration or on landing. If all of the 
body’s potential and kinetic energy were converted into thermal 
energy within the body itself, it would be more than sufficient to 
vaporize the entire body, unless the body were entirely con- 
The survival of 
meteorites, however, is an obvious indication that not all of the 
energy goes into the body. 

Actually, the body’s initial energy is transformed, through the 
mechanism of gasdynamic drag, into thermal energy in the air 
around the body, and only part of this energy is transferred to the 
body as heat. The fraction of the original energy that ap- 
pears as heat in the body depends upon the characteristics of the 
flow around the body. 

Effects During Penetration. Before entry into a 
planetary atmosphere, a body’s motion is governed only by the 
principles of celestial mechanics; 


structed of some material such as carbon. 


Gasdynamic 


i.e., its motion is governed 
only by its own inertia and by gravitational forces. However, 
such motion is modified when the vehicle enters the outer fringes 
of a planetary atmosphere. In the rarefied upper atmosphere 
where the molecular mean free path is large compared with the 
body size, the flow is of the so-called free-molecular type and a 
gasdynamic drag force begins to act in a direction opposite to the 
vehicle’s motion. In the free-molecule flow region, appreciable 
lift forces (i.e., forces normal to the direction of motion) are not 
attained (at velocities of the orbital or escape magnitude). In 
the denser lower atmosphere, where the flow is of the continuum 
type, both lift and drag forces affect the vehicle’s motion. 

The effects accompanying atmospheric entry that are mani- 
fested in the vehicle are a reduction in velocity accompanied by 
appreciable deceleration loads and by appreciable heating rates. 
The deceleration, and consequent velocity reduction, is due to 
aerodynamic drag, which varies with air density (p) and the 
square of the vehicle velocity (u?): 


Deceleration ~ pu? 


Now, as the vehicle approaches the planet it first encounters an 
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atmosphere of very low density. As it penetrates further, the 
density increases rapidly; and, due to drag, the velocity begins 
Thus the deceleration is the product of two quan- 
Initially, the 


however, at some point, the velocity be- 


to decrease 
tities, one increasing and the other decreasing 
deceleration increases; 
gins to decrease more rapidly than the density increases, resulting 
in a maximum deceleration with subsequent decreasing decelera- 
tion in most cases. This is illustrated in Fig. 2. For entry from 
space, these decelerations may be very appreciable; for example, 
entry at right angles into the earth’s atmosphere at escape velocity 
involves a maximum deceleration of over 300 g. 


ALTITUDE 
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VELOCITY, uv DECELERATION HEATING RATE 


Fig. 2 Changes during atmospheric entry 


Gasdynamic heating of the vehicle surface during entry varies 
in a somewhat similar fashion. The heating rate varies approxi- 
mately with the atmospheric density and the cube of the vehicle 


veloc ity 
Heating rate ~ pu? 


Here, also, the increasing density and decreasing velocity during 
entry interact to produce a maximum heating rate. This is also 
illustrated in Fig. 2 

These « xpressions indicate that both deceleration and heating 
are most severe when there is a combination of high atmospheric 
density and high vehicle velocity, 3.C., when high veloc ities are 
allowed to persist to low altitudes. 
in which 


Factors Affecting Penetraiion Characteristics. The way 


velocity changes with altitude during entry, the deceleration 
forces, and the surface heating rates depend upon several factors: 
Atmospheric approach (Fig. 3 
Velocity of approach, u, 
Angle of approach, @ 
o Vehicle characteristics (Fig 
W : 
Mass-drag ratio, (W = 
Cp)A 


drag coefficient, 


| frontal area 
(Atmospheric characteristics ( Fig. 5) 
Surface density, ps: 


Density variation with altitude, p = f(h 


The velocity and angle of atmospheric approach, illustrated in 
Fig. 3, are primary factors in determining the severity of the entry 


deceleration and heating. As indicated previously, combina- 


tions of high velocity and low altitude produce the most severe 


entry conditions. These are most apt to occur when the initial 
velocity is high or when the entry is perpendicular (6 = 90 deg). 
On the 


angle tends to restrict higher velocities to higher altitudes 


other or a shallow entry 


The 


initial entry velocity is determined by the planet’s gravitational 


hand, a lower initial velocity 


characteristics and by the type of entry, i.e., from satellite orbit, 
directly from space, ete. Since usually one is not able to choose 
the initial entry velocity (except by rocket braking), adjustment 
Kf- 


fectively, high velocities can be restricted to the upper atmosphere 


of the angle is one means of reducing the severity of entry. 


either by a lowered initial velocity or by a shallow entry-angle 30 
that deceleration takes place high in the atmosphere. 
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Fig. 3 Entry into a planetary atmosphere 


Similarly, deceleration high in the atmosphere can be affected 
by the use of a body having a small mass-drag ratio and/or some 
aerodynamic lift. A reduced mass per unit drag allows decelera- 
tion in the less dense region of the atmosphere, while aerodynamic 
lift decreases the rate of descent and, thus, allows time for de- 
For example, the slender 
body shown in Fig. 4(a) has a high mass-drag ratio and would 
penetrate (at high velocity) further into the atmosphere than the 


celeration high in the atmosphere. 


blunt body with a low mass-drag ratio illustrated in Fig. 4(b 
If the latter body is oriented in the position shown in Fig. 4(c), a 
lift foree would be developed which would force it into a more 
shallow path and hence restrict deceleration to a high altitude 


0. SLENDER,LOW-DRAG BODY 


ee... 


o 
Biro to } 
b BLUNT 
a 
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i 


BODY 


c. LIFTING BODY 


Fig. 4 Aerodynamic forces on various bodies 


The physical and chemical characteristics of the planetary 
atmosphere also determine the entry characteristics. For a first 
approximation to the gasdynamic forces and heating, a knowledge 
of the density variation in the atmosphere is sufficient. For hy- 
drostatic equilibrium in an isothermal atmosphere, it can be 
shown that the atmospheric density at altitude, in ratio to the 


sea-level density, is given by 


where e is the base of natural logarithms, h is altitude, and a@ is a 
coefficient dependent upon the temperature and composition of 
the atmosphere. This equation indicates that the density de- 
While the 


planetary atmospheres are far from isothermal, an expression of 


creases exponentially as the altitude increases 


this type represents the earth’s atmosphere reasonably well. 
The approximate density variation in the earth’s atmosphere and 
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Fig. 5 Density variation in planetary atmospheres 


the estimated density variation in the atmospheres of Venus and 
Mars are shown in Fig. 5 [1].!_ The atmosphere of Venus, esti- 
mated to consist of about 10 per cent nitrogen and 90 per cent 
carbon dioxide, is somewhat earth’s at- 
mosphere but varies in a similar way with altitude. The at- 
mosphere of Mars, estimated as about 95 per cent nitrogen and 


more dense than tie 


5 per cent carbon dioxide, is appreciably less dense at sea level 
but varies much more gradually with altitude and is actually 
more dense at high altitudes. The more gradual density varia- 
tion in the Martian atmosphere effectively makes it “softer’’ and 
consequently reduces the severity of entry. 


Dynamics of Atmospheric Penetration 


As illustrated in Fig. 1, there are three general types of at- 
mospheric penetration: 


1 Relatively steep descent path, characteristic of direct entry 
from space, “dump” 
ballistic rocket. 

2 More gradual descent path, characteristic of a satellite orbit 
decay. 


from a satellite orbit, and re-entry of a 


3 Very gradual descent, characteristic of a boost-glide rocket 
or a lifting descent from a satellite orbit. 


While the equations of motion for atmospheric penetration re- 
quire a numerical or machine computation for exact solution, ap- 
proximate analytic solutions are possible for the three cases listed 
above and have been described elsewhere [1] 

Direct Entry. For the case of a relatively steep (@ > 5 deg) direct 
entry on a linear path, a simple approximate solution may be ob- 
tained for a nonlifting vehicle. If the gravitational force is 
neglected in comparison to the decelerating aerodynamic force, 
if the atmosphere is assumed isothermal, and if the drag co- 
efficient is assumed constant, the velocity variation is given by 
(2, 3, 4, 5] 

CoA Psi 


Wsin@ 2a 


Comparison of this approximate equation with the results of exact 
machine computations indicates good agreement for cases where 


! Numbers in brackets designate References at end of paper 
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the entry angle is greater than about 5 deg and in the velocity 
range where deceleration loads and gasdynamic heating are im- 
portant considerations. (For small entry angles, Cases 2 and 3 
above, other approximations apply.) After the major decelera- 
tion has occurred, the vehicle’s path bends over into a vertical 
descent; and the velocity is essentially a local terminal-fall 
velocity. 
Fig. 6 illustrates the velocity variation with altitude for direct 
penetration of the earth’s atmosphere (a = 4.15 & 107° per ft) 
CpA, 
W sin 0 


for several values of the drag-mass-angle parameter, 


400 


ALTITUDE, KFT 











Fig. 6 Velocity variation during direct penetration of the earth's 
atmosphere 


One interesting facet of these results is that the major decelera- 
tion of the body takes place in a stratum of the atmosphere of 
constant thickness. 


The level of this stratum in the ~*mosphere 


is, however, dependent upon the value of the drag-:mass-angle 


parameter. Large values of this parameter cause deceleration 
high in the atmosphere, while smaller values delay decelera- 
tion until the lower atmosphere is reached. Thus the use of a high- 
drag, light body at a small entry angle enables deceleration at 
high altitudes. Since severe deceleration loads and gasdynamic 
heating rates occur when high velocities persist down to low 
altitudes, the use of a light, high-drag body at a shallow entry 
appears advantageous. 

The deceleration during a relatively direct atmospheric pene- 
tration is found to rise to a maximum value and then decrease 
during descent (see Fig. 2). The maximum value of the de- 
celeration occurs when the velocity has been reduced to about 


61 per cent of its initial value. 
= (0.607 


The maximum deceleration is 


( du au,? sin @ 
dt max . Ze 


It is interesting to note that the maximum deceleration during 
direct entry is independent of the drag-mass characteristics of the 
body. It is dependent only upon the path angle, the initial ve- 
locity, and the atmospheric characteristics. For a given entry 
velocity, entry angle, and atmosphere, the maximum deceleration 
is independent of the body characteristics; however, the altitude 
of maximum deceleration is dependent upon the mass-drag charac- 
teristics of the body. 
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Fig. 7 Velocity and deceleration during direct entry into the earth's 
atmosphere from space at various angles W/(CpA.) = 10 Ib/sq ft 


Fig. 7 illustrates the effect of entry angle on the variation of 
velocity and deceleration during penetration of the 
atmosphere from space by a body having the mass-drag charac- 
teristics of the Vanguard scientific satellite (W = 21.5 lb, diam- 


earth’s 


eter 20 in., ~ 10 lb/sq ft). The smaller entry angle is 
CpA, 
deceleration higher in the atmosphere and a 


seen to result in 


smaller peak deceleration. If a lighter body were used at these 
entry angles, these curves would be shifted to higher altitudes 


A heavier bod) wuld result in a shift to lower altitudes 
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Fig. 8 Velocity and deceleration during direct entry from space at 
@ = 90 into three planetary atmospheres X/(CpA.)/ = 10 Ib/sq ft 


rhe variation of velocity and deceleration of the same body is 
shown in Fig. 8 for vertical penetration of the atmospheres of 
Earth, and Mars [1]. The different gravitational at- 
traction accounts for the different initial velocities. The effect of 
atmospheric density variation is apparent in the shape and posi- 
tion of the curves 


Venus, 


Che more dense atmosphere of Venus results 
in deceleration at a higher altitude; however, the velocity varia- 
tion with altitude and the maximum deceleration is about the 
same as for Earth because of the similar density variation (see Fig 
5). However, the more gradual variation of density in the 
Martian atmosphere results in a more gradual variation of veloc- 
ity with altitude and a lower peak deceleration. 

Orbit Decay. In the case of the more gradual descent path 
characteristic of the decay of a satellite orbit (see Fig. 1), the as- 
sumptions used in the previous section are no longer applicable. 
The path is far from linear, and the path angle changes during 
descent. 
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In the region of appreciable deceleration and heating, 
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the path angle is of the order of 3-4 deg. An approximate 
analytic solution for this case [1] yields results somewhat similar 
to those for direct entry in that the altitude of deceleration in the 
atmosphere is dependent upon the mass-drag characteristics of 
the body and the maximum deceleration is independent of these 
characteristics. Velocity variation with altitude is shown in Fig. 


9, with 


’ as parameter. 
Dé*e 

this parameter since it cannot be specified and since it varies dur- 
ing descent—it is an output from rather than an input to the 


The deceleration varies in a somewhat similar manner 


(The path angle does not appear in 


problem. ) 
to that in the previous case, rising to a maximum and then falling 
off. Here the maximum deceleration occurs when the velocity 
has been reduced to 


and has a relatively low value because of the low path angle. For 
example, the decay from a circular satellite orbit into the earth’s 
atmosphere involves a maximum deceleration of about 9g. 
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Fig. 9 Velocity variation during penetration of the earth's atmosphere 
by orbital decay 


It should, perhaps, be emphasized that the orbital decay de- 
scribed here refers to the last stage of decay where heating and 
about the last 2000 miles and last 
few minutes of the satellite lifetime. 


deceleration are appreciable 
This phase is preceded by 
a much longer period, many revolutions, in which the satellite 
orbit gradually decays. This initial phase of decay consists of an 
extremely gradual spiral during which the orbit becomes more 
circular. The rate of energy loss by the vehicle due to aerody- 
namic drag is small enough so that the vehicle’s kinetic energy 
and potential energy adjust themselves to a momentary “ equi- 
librium’’ orbit. In this process, the potential energy decreases 
and the kinetic energy increases—with a net decrease in total 
energy. Thus the satellite’s velocity actually increases in the 
initial phases of orbital decay. 

Lifting Descent. A lifting descent involves a still more gradual 
atmospheric penetration (see Fig. 1). Here again, the path 
angle adjusts itself to the forces acting on the vehicle and turns 
out to be quite small—of the order of a few tenths of a degree. 
The velocity variation with altitude, as indicated by an ap- 
proximate analytic solution [1], is shown in Fig. 10 as a function 
CA. LCpA, 

W° Dw 
drag-mass parameter by the lift-drag ratio, L/D. 


of a lift-mass parameter, , which is related to the 


It is seen that 
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here again the level of deceleration in the atmosphere is deter- 
mined only by the mass and aerodynamic-force characteristics of 
the body. In this case, however, the deceleration does not go 
through a maximum value, but asymptotically approaches a 


value 
( du g 
Gsun L/D 


Since aerodynamic shapes with lift-drag ratios as high as five 
appear feasible, it is clear that relatively small decelerations can 
be obtained in a lifting descent. 
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Fig. 10 Velocity variation during penetration of the earth's atmosphere 
by a lifting vehicle 


Comparisons. In comparing the three types of atmospheric 
penetration just discussed, it is necessary to keep in mind the 
fact that the more gradual descents involve longer times and 
greater ranges than a steeper descent (see Fig. 1). For example, 
descent through the ‘‘Corridor of Continuous Flight’’ [6] can be 
accomplished by the three types of penetration described with 
essentially the same velocity-altitude variation, as shown in Fig. 
11. However, the direct descent would involve a distance of 
The 
orbital decay might involve a range of about 2000 miles and 5 or 
10 minutes. The lifting descent would extend over 5000 to 10,- 
000 miles and about two hours. 


only up to a few hundred miles and about one half minute. 


A gradual descent involves 
velocity reduction and consequent energy dissipation over a long 
period of time. Thus deceleration loads and aerodynamic heat- 
ing rates are reduced by lengthening the time of deceleration. The 
maximum deceleration loads experienced during penetration of 
the earth’s atmosphere, with initial entry at satellite velocity, is 
shown in Fig. 12 for direct descents at various entry angles, for 
orbital decay, and for descent with various degrees of aerody- 
namic lift. The maximum deceleration loads for direct entry from 
space at escape velocity would be twice as high as those shown for 
direct entries at satellite velocity. It is clear that a steep direct 
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Fig. 11 Velocity variation with altitude for three types of atmospheric 
penetration 
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Fig. 12 Maximum deceleration into the earth's atmosphere at orbital 
velocity 


entry at either velocity would probably not be humanly tolerable 
and that the atmospheric penetration by a manned vehicle would 
require a relatively shallow descent path. The descent paths for 
vehicles with aerodynamic lift involve quite tolerable deceleration 
less than one g. The maximum deceleration for various 
types of entry into the atmospheres of Venus, Earth, and Mars 


are given in Table 1. 


loads 


Heating During Atmospheric Penetration 
As described earlier, the reduction of the vehicle’s kinetic and 
potential energy during descent is accompanied by an increase in 


Maximum deceleration experienced during various types of atmospheric penetration 


(Values given in Earth g's) 


Direct entry at escape 
Planet velocity 
6=5 20 90 
Venus 28.6 
Earth 28.3 324 
Mars 1.6 : 
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Direct entry at orbital 
velocity 


Entry by 
decay from 
satellite 
orbit 


Entry of lifting vehicle 
at orbital velocity 


20 90 : 5 


163 8.9 8 18 
162 9.5 2 
9.2 1.4 076 
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thermal energy in the surrounding air, some of which is transferred 
to the vehicle’s surface. Now the fraction of the vehicle’s lost 
kinetic and potential energy which reaches the vehicle surface as 
heat is of primary concern to the designer. This fraction, or 
conversion efficiency, depends upon the vehicle's shape and on 
its velocity and altitude—and ultimately on the mechanism of 
heat transfer between tie hot gas and the vehicle surface. At 
very high altitudes, in the free-molecule flow regime, where the 
heat energy is developed directly at the vehicle’s surface, as 
much as one half of the vehicle’s lost energy can appear as heat 
in the body. At lower altitudes, in the continuum flow regime, 
thermal energy appears in the air between the shock wave and 
the body. Heat is transferred from this hot air to the body by 
conduction and convection through a viscous boundary layer; 
radiation from the hot gas may also contribute appreciably to the 
surface heating. 


The rate of change of the vehicle’s kinetic energy per unit mass 


WW 2 


is 


pu 
The portion of this reaching the body as heat, i.e., the heating 
rate per unit Mass, can be expressed as, approximats ly, 


q l 


yr Cuan pu 


where Cy, is a dimensionless heat-transfer coefficient (Stanton 


number) based on frontal area A,, and 7,» is the recovery factor 
tatioing these two equations gives the instantaneous “‘efficiency’’ 
of converting the vehicle's kinetic energy to heat energy in the 
body 

Cun 


f= 
. ( 


D 
The value of this conversion fraction f de pends upon the heat- 
transfer mechanism, the velocity, the altitude, and the vehicle’s 
geometry. A typical variation of f for a relatively blunt body at 
a velocity of 35,000 ft/sec at various altitudes is shown in Fig. 13 
It is seen to vary from '/, at very high altitudes down to about 


0.01 at 100,000 ft At still lower altitudes it increases again due 


400 ; 


LE HEATING ~ 





FRACTION ENER 
Fig. 13 Conversion efficiency (kinetic to thermal energy) as a function 
of altitude 
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to a turbulent boundary-layer condition and to appreciable radia- 
tion from the hot gas. At a lower velocity, e.g., satellite velocity, 
the radiation contribution becomes small and can usually be 
neglected. 

The maximum heating rate, for the several types of penetration 
described in the previous section, occurs when the velocity has 
The 
peak heating rate, therefore, occurs before the point of maximum 


been reduced to about 85 per cent of the initial velocity. 


deceleration load. 

The surface heating rate is dependent both upon the rate of re- 
duction of the vehicle’s kinetic energy and upon the conversion 
At high altitudes where, although f is high, the rate 
of energy reduction is low, the heating rates will be low. Similarly, 


efficiency. 


at moderate altitudes, where f is low, the rate of energy reduction 
is high and high surface heating rates will occur. This is il- 


lustrated in Fig. 14. 
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Fig. 14 Variation of rate of reduction of vehicle's kinetic energy and 
surface heating rate 


The total heat input to the body during atmospheric penetra- 
tion is obtained by integrating the heating rate over the time of 
descent. If the conversion factor f were constant, the total heat 
load would be independent of the type of entry, since a given frac- 
tion of the body’s initial energy would be converted into heat in 
the body. However, the actual variation of f (Fig. 13) allows one 
to minimize the total heating by positioning the maximum rate of 
energy reduction at the minimum value of f. In the case of the 
examples illustrated in Fig. 14, the total heat load for the vehicle 
decelerating in the upper atmosphere would be greater than that 
for the vehicle decelerating at a lower altitude. The former 
would involve lower heating rates over a longer period of time, 
but with a greater accumulated heat load. 

If the heating rate is sufficiently low, it can be rejected by ra- 
diation from the heated surface. This type of cooling is, of 
course, limited by the temperature limitations of the surface 
material. When all of the incoming aerodynamic heating is re- 
jected by thermal radiation from the surface, we may write a 
simple heat balance 


q a 
~ eT, 


which says that the incoming heat from hot air surrounding the 
body is proportional to the surface emissivity and the fourth 
power of the surface temperature. This temperature is the so- 
called radiation equilibrium temperature, and serves as a con- 
In Table 2 
the maximum heating rate, in terms of this radiation equi- 


venient parameter for the expression of heating rates. 


librium temperature, is given for several types of penetration into 
the atmospheres of Venus, Earth, and Mars for a body with a 
heated surface area twice its frontal area. In some cases, it is 
apparent that the temperatures are within the range of current 
materials. In other cases, the heating rates correspond to tem- 
peratures above the range of current materials. 

Whether it is desired to minimize the peak heating rate or to 


minimize the total heat load depends upon the penetration 
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Table 2 Maximum heating rate experienced during various types of atmospheric penetration for 


CyA 


(Values expressed as equivalent radiation equilibrium temperature, deg R) 


Direct entry at escape 

Planet velocity 
@=5 20° 90° @ = 5° 
4700 
5000 
2400 


Venus 
Earth 
Mars 


5500 
5900 
2800 


6300 
6800 
3200 


3600 
3900 
1900 


Direct entry at orbital 
velocity 


4200 
4600 
2200 


Entry by 
decay from 
satellite 
orbit 


Entry of lifting vehicle 
at orbital velocity 


20° 90° ~ = } 2 5 


D i 
2300 
2500 
1500 


4900 
5200 
2500 


3400 
3600 
1900 


2100 
2300 
1350 


1900 
2000 
1200 


* Note that the maximum heating rate, unlike the maximum deceleration, is dependent on the mass- 


drag parameter 


technique chosen. If one chooses to reject heat by thermal radia- 
tion, it is desirable to minimize heating rates and hence surface 
temperatures. On the other hand, if heat is being absorbed by a 
vaporizing surface, it is desirable to minimize the total heat load. 
These considerations are discussed more completely in the next 


section. 


Penetration Techniques 

While a great variety of techniques for successful atmospheric 
penetration have been suggested, a synthesis of the principles 
outlined in the preceding sections suggests two general classes: 


1 The first involves the use of a vehicle having a relatively 
small mass-drag or mass-lift parameter, or a shallow entry path, 
so as to effect deceleration high in the atmosphere. In this case, 
the aerodynamic heating rates are low enough so that the heat 
may be rejected by thermal radiation from the vehicle surface. 
This corresponds to a 
tinuous Flight [6] 


descent through the Corridor of Con- 
since the descent is slow enough so that 
steady-state conditions apply. Examplk s of this type of descent 
are the descent of a glide rocket, the orbital decay of a Vanguard 
satellite [7], or the direct entry from space of a relatively light, 
high-drag body. These three example s involve successively 
greater deceleration loads—thus a manned descent would proba- 
bly require either a lifting vehicle or a very shallow descent. 

2 The second general class involves the use of a vehicle whose 
mass-drag characteristics and entry path are adjusted so as to 
place the point of maximum energy reduction at the point where 
the efficiency of conversion of vehicle energy to surface heat is a 
minimum, i.e 


, so the total heat load is minimized. In this case, 


the heat is absorbed by the body surface. An example is the direct 
descent of an unmanned vehicle from space [8]. The technique 
of heat absorption may consist of a solid heat sink, internal cool- 
ing, or a vaporizing surface. This type of atmospheric penetra- 
tion generally involves higher deceler ition loads and is most ap- 


plicable to unmanned vehicles. 


While it does not appear possible to specify the most desirable 
system of atmospheric penetration for a given mission without a 
detailed design study, one can indicate general areas of feasibility. 
For example, Fig. 15 indicates a estimate of the 
penetration systems applicable to various situations. 


schematic 
Here are 
plotted the entry angle and the mass-drag ratio. Small values of 
these variables yield relatively low heating rates and radiation 
cooling is feasible. Large entry angles and large mass-drag ratios 
mean high heating rates (but lower total heat loads), and heat 
absorption schemes are applicable. Thus the regions of radia- 
tion cooling and heat absorption are indicated. Deceleration 
loads are primarily dependent upon entry angle (for a given entry 
velocity). Manned descent is therefore restricted to the smaller 
angles. It will be noted that both radiation and heat absorption 
techniques are applicable to both manned and unmanned at- 


mospheric penetration. 
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Fig. 15 Approximate regions for various types of entry from space 
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Fig. 16 Cooling-system weight requirements for direct entry from space 
at escape velocity 


One method of deciding which type of cooling scheme to use is 
the examination of the weight which must be charged to the 
solution of the heating problem. This weight may be either 
that of a heat absorbing surface or the weight of drag brakes or 
lifting surfaces used to bring about the reduced heating rates 
characteristic of high-altitude deceleration. Such cooling system 
weight, in ratio to total weight, is shown in Fig. 16 as a function 
of the mass-drag parameter. At low values of the mass-drag 
parameter (large, light bodies) weight is required for drag brakes 
or lifting surfaces and surface heating is handled by radiation. 
At large values of the mass-drag parameter (small, heavy bodies) 
weight is required for heat absorption. The shape of the curve is 
seen to be somewhat similar to the conversion efficiency (Fig. 
13)—although here two minimums appear. One minimum weight 
corresponds to the point where f is a minimum, i.e., where the 
The 


other minimum corresponds to the point where the mass-drag 


sum of the convective and radiative heating is a minimum. 


parameter is just small enough so that radiative cooling is possible 
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and yet not large enough to require a large weight in drag brakes 
or lifting surfaces 
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Thermodynamic Properties of Helium at 
Low Temperatures and High Pressures 


The thermodynamic properties of helium have been compiled and correlated for a tem 
perature range from 3.0 to 20° K for pressures to 100 atmospheres and for specific 
volumes from 5 to 800 liters per kilogram. 
temperature-entropy and the enthalpy-entropy co-ordinate systems and include pressure, 


The properties are presented on both the 


temperature, volume, entropy, and enthalpy 


Introduction 


= INTEREST in compressed liquid helium has 

resulted in considerable new experimental data of the thermody- 
namic properties of helium. Investigations of liquefaction or 
refrigeration cycles and the design of heat exchange equipment 
indicate the need for a correlation of the thermodynamic proper- 
ties and their presentation in an accurate and convenient manner. 
The most comprehensive investigations of the thermodynamic 
properties of normal helium I and II were published over a 
period of the last thirty years by W. H. Keesom, and a culmina- 
tion of most of this work was published by Keesom [1]? in 1942. 
Since this time, numerous investigators have determined con- 
siderable additional data, particularly values for compressed 
liquid in the temperature range below 20 K, an area not covered 
extensively by Keesom. These recent investigations in this low 
temperature range have been published as separate and some- 
times isolated articles of very narrow interest and in some cases 
The authors have made 
a compilation and a critical correlation of the present data on the 


the results are not in readily usable form. 


thermodynamic properties of helium below 20 K and have pre- 
sented them in the form of two convenient charts. The resulting 
diagrams include the properties of pressure, volume, tempera- 
ture, enthalpy, and entropy for pressures up to 100 atmospheres, 
for specific volumes from 5 to 800 liters per kilogram, and for tem- 
peratures between 3.0 and 20 K. The first of these charts uses 
temperature-entropy co-ordinates which will be generally useful 
in analysis of thermodynamic systems and processes. The second 
chart presents the data with enthalpy-entropy co-ordinates. 
This second chart was found necessary since the slopes of the 
volume, pressure, and enthalpy lines approach a common value 
in the compressed liquid region on the temperature-entropy plot, 
making interpolation difficult and inaccurate. This graph is not 
extended to low pressures and high temperatures as the tempera- 
ture-entropy diagram published by J. Zelmanov [2] is adequate 
for these regions. 


Data 


The greatest quantity of information prior to 1944 on the ther- 
modynamic properties of helium in the region being considered 
comes from two sources. The first of these is a compilation by W. 
H. Keesom [1] of the theoretical and expt rimental investigations 
performed at the University of Leiden, Holland. This work is 

! Results of a study made under contract with the Bureau of Aero- 
nautics, Department of the Navy at the National Bureau of Stand- 
ards, Boulder Laboratories. 

? Numbers in brackets designate References at end of paper 

Contributed by the Heat Transfer Division of THe AMERICAN 
Society oF MECHANICAL ENGINEERS and presented at the ASME- 
AIChE Heat Transfer Conference, Storrs, Conn., August 9-12, 1959. 

Nore: Statements and opinions advanced in papers are to be 
understood as individual expressions of their authors and not those of 
the Society. Manuscript received at ASME Headquarters, May 20, 
1959. Paper No. 59-—-HT-19 
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quite extensive and covers the temperature range from 0 to 400 K 
and pressures from 0.001 to 200 atmospheres. Secondly, an in- 
vestigation of the Joule-Thomson effect in helium by J. Zelmanov 
2, 3] adds considerable data at low temperatures. This work was 
undertaken because of discrepancies found in Keesom’s values at 
temperatures below 20 K, and at pressures between 0.5 to 60 at- 
mospheres. 

More recent data have been presented by R. W. Hill [4, 5] in 
the region from the critical pressure of 2.2 atmospheres to 100 at- 
mospheres and for temperatures from 2.5 to 20 K. The data 
given are values of temperature, specific volume, and entropy. 
Data from O. V. Lounasmaa [6], which only recently have be- 
come available, list experimental values of pressure for the data 
presented earlier by Hill. Several articles on the temperature 
scale add to the data on vapor pressures of helium of which some 


use has been made. A summation and extension of the experi- 


mental data of the solid for a pressure range of 25 to 225 atmos- 
pheres and for temperatures from 0 to 6 K, have been published 
by C. Domb and J. 8. Dugdale [7]. In addition, Berman and 
Mate [8] and Van Dijk and Durieux [9] have presented extensive 


data on temperature, vapor pressure, saturation densities, and 
heats of vaporization. 


CONSTRUCTING THE CHART 
Reference Values for Entropy and Enthalpy 


It was first necessary to establish a fixed point for values of 
entropy and then to refer the data of the several investigators to 
this fixed point In addition, the temperature scales used by the 
earlier investigators were corrected to agree with a presently ac 
ceptable temperature scale designed as 7';£, proposed by Clem- 
ent, Logan, and Gaffney [10]. The normal boiling point of 
liquid helium was selected as the fixed reference point, that is, 
4.21 K and one atmosphere pressure. The value of entropy at 


this reference point was calculated by the relation, 


eT =4,.21°K 
ove fe 
ag a 
J T=0°K 1 


where S = Ocal gm deg K at 7 =OK. Then, according to the 
Nernst-Simon statement of the Third Law of Thermodynamics 
[1, 11], the entropy of the solid at 0 K and 25 atmospheres pres- 
sure and the saturated liquid for vaporization at 0 K and zero 
pressure are identical and equal to zero. Thus the values of en- 
tropy for the saturated solid from Domb and Dugdale [7], which 
were reported with 0 K and 25 atmospheres at zero entropy as 
the reference point, could be compared with the values of entropy 
from the other references. The values of C,.: from the data of 
Kramers, Wasscher, and Gorter [12] were used to calculate the 
entropy differences below 1.9 K, while the values of Cy. between 
The 
entropy difference between saturated liquid in equilibrium with 
its vapor at 0 K and at 4.21 K was determined as 0.830 cal/gm K, 
which is in agreement with the value reported by Lounasmaa 
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1.8 K and 5.05 K are from data of Hill and Lounasmaa. 
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Phe fixed point tor enthalpy was also selected at the normal 


his value of 2.39 cal/gm was 


. oll id 
fom +f" var, 
7 / 


° 
using C,.. values as discussed in determination of the 


boiling point of liquid helium 


determined by the relation, A h = 


satu 


ration entropies. Values of vapor pressure were from Clement, 


ind Gaffney 


from Berman 


Logan ind values of specific volume at saturation 


Construction of Isobars 
Lines of constant pressure were ¢ ileulated for the charts from 
oS /ol oP /ol 


Only the pressure differences ilong the isometrics of 5. 5.5, 6, 7 


Hill’s data using the Maxwell relation, 


1 


and 8 liters per kilogram were used for this purpose as the Zel- 


manov data were found to adequ itely describe the isometrics for 
larger values up to 800 liters per kilogram. The calculated 
values for pressure were also compared with values from Lounas- 
maa when they became available Differences between these cal- 
culated values and Lounasmaa’s experimental values were less 
than 5 per cent, confirming the compatibility of these data. Since 
the calculated 


values of pressure of the authors are considered in error when in 
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the calculation of slopes is inherently inaccurate, 
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Thermodynamic properties of helium. 


conflict with the experimental pressure data from Lounasmaa 
In addition to the Lounasmaa data, values for pressures from 0.5 


to 2 atmospheres are from Zelmanov. 


Construction of Isenthalpic Lines 


Values of enthalpy at saturation were determined as discussed 


} . 
Love 


Starting from these saturation values, the values of en- 
thalpy in the compressed liquid egion were then found using the 
relation, dh = 7 ds + 


the normal boiling point, using small increments of pressure and 


vy dP, first along an isentropic path from 
average values of specific volume. Having established enthalpy 
values along this isentropic line at pressures to 100 atm the 
other values of enthalpy were determined along the constant 
pressure paths When the Lounasmaa data became available, 
these calculations were again checked using the new values for 
pressure. Below 2 atmospheres pressure, the data of Zelmanoy 


were used in addition to the Lounasmaa data 


Values for the Isomeric Lines 


Specific volumes from 5 to 70 liters per kilogram were from 
Lounasmaa and Hill and from 70 to 800 liters per kilogram from 
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. exe . = a - 
Zelmanov The values for specific volume 7 to 70 liters per 


kilogram from Zelmanov, are also in excellent agreement with 


Hill. 


existed the data from Lounasmaa and Hill wer 


those from Lounasmaa and However, where differences 


used 


Properties of the Saturated Vapor 


The heats of vaporization from Berman were used to deter- 
and saturated 


dP 
dT 


mine the entropy differences of saturated vapor 


liquid. By using the Clausius-Clapeyron equation, 
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h 


» data from Berman were also used to determine 


T(v, . v; 
densities of the gases at saturation. 
Properties at Solid-Liquid Equilibria 


The Domb and Dugdale data were used for these boundaries. 


Lambda Transition 


Properties along this boundary were taken from Keesom 
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Conclusions 


rhe authors have carefully compared all of the data available 
in the region investigated, making extensive calculations to in- 
sure the compatibility of the data from the various sources. In 
all but 


various sources differed less than one per cent 


1 limited area concurrent values determined from the 


This limited area 
of larger discrepancies was between 2 and 60 atmospheres pres- 
sure for temperatures from 12 to 20 K. In this area many of the 
values from Zelmanov [3, 4] and Lounasmaa [6] differed more 
than one per cent, ranging to as much as three per cent at the 
higher temperatures. These discrepancies were resolved pri- 
marily by comparison of these data to values in adjacent areas. 
As a result of a careful review of the experimental apparatus and 
the techniques used, and making use of each investigator's 
estimate of the probable accuracy of his own data, it is believed 
that the values presented on the charts are accurate to within 5 


per cent and in most cases to within 3 per cent. 
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The one-dimensional steady-state solution for radiation from one side of an annular 
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Introduction 


I. HAS LONG BEEN KNOWN that the rate of heat trans- 
fer from a hot body to a cooler fluid may be increased by extend- 
ing the body surface. Furthermore, the addition of fins extending 
into the fluid will cause an increased heat flow not directly pro- 
portional to the increase in body surface, since the temperature 
near the tip of the fin approaches the fluid temperature and 
causes a heat rate less than that which might be anticipated. 

All of the previously published literature considering heat trans- 
fer from an extended surface has considered only the convective 
type of heat transfer. Convection from straight fins was con- 
sidered by Parsons and Harper “{1],?, Harper and Brown [2], 
Avrami and Little [3], and Carrier and Anderson [4]. Convec- 
tion from annular fins was the subject of papers by Harper and 
Brown [2], Carrier and Anderson [4], Murray [5], Gardner [6], 
and Keller and Somers [7 }. 

Only recently, radiation heat transfer from fins has become of 
practical engineering importance; especially in the field of 

1 This paper is based on an MS thesis submitted by R. L. Chambers 
t’ the Mechanical Engineering Department of the University of 
Fittsburgh on the Westinghouse Advanced Mechanics Fellowship 
Program. 
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fin has been computed for values of the two design parameters 0) 


2.0 and 1.001 < p= 


eo@ 3 
K (r, — 7) 
r , - | . . . 

< 15.0. The solution supplies design inforn-ation needed 


‘ 
for satellite thermal-power dissipating surfaces 


auxiliary power units for satellite vehicles or high altitude mis- 
siles. At high altitudes such vehicles dissipate heat only by radia- 
tion, since atmospheric convection is impossible. In the design of 
such power units, the principle design factor of weight per unit 
power is controlled by the heat dissipation rate at the heat-rejec- 
tion surface of the engine, i.e., by the temperature distribution of 
the surface, since the heat rejection rate is controlled by the 
fourth power of the absolute temperature [8], i.e., the Stefan- 
Boltzmann law. Weight reductions are accomplished by finning 
techniques, and thus the dissipation of heat by radiation from 
fins has become important. To date, estimates of the fin ef- 
ficiency for these applications have been made by linearizing 
the radiation rate, i.e., by replacing the fourth-power law by an 
equivalent convection constant multiplied by the first power of 
the absolute temperature, and then using the available data al- 
ready cited [6]. That this procedure will introduce significant 
errors in the computed heat rate, as great as 60 per cent, will be 
shown by a comparison of the results of such a procedure with the 
results of this paper. The work reported here discusses the solu- 
tion of a typical fin geometry applicable to certain auxiliary power 
units and cooled only by radiation from the fin surface to space. 


Statement of the Problem 


The problem is to obtain the steady-state, one-dimensional 
temperature distribution and fin efficiency for a flat mushroom 
The fin 
conduction 


annular fin with variable size and material properties 
geometry is shown in Fig. 1. The heat enters by 


through the peg and leaves by radiation from the exposed surface 





Nomenclature 


= two-side heat-transfer coefficient 

for convection, Btu/hr ft? deg 

F Btu/hr 

one-side heat-transfer coefficient 

for convection, Btu/hr ft? deg F 

thermal conductivity of fin ma- 
terial, Btu/hr ft deg F 


: it 
total number of concentric ele- 


dimensionless 


ments that subdivide the fin 
) ae ’; 


|2h lh’ f To — 5 
= a 
Vis = Vis 
n an integer 
Q heat flow rate, Btu/hr 
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heat flow rate of a fin at the con- 
stant temperature of the base, 


radius at any point in the fin, ft 
inside radius of the fin, ft 
outside radius of the fin, ft 
width of a concentric fin element, 


emissivity, dimensionless 
Stefan-Boltzmann 

0.1714 K 10-* Btu/hr ft? deg R‘4 
fin thickness, ft 


fin efficiency, dimensionless 

dimensionless ratio of r,/r; 

temperature at some point in the 
fin, deg R 

temperature at the fin base or r,, 
deg R 

dimensionless ratio of 


ratio of radii = 


A'(y 
dimensionless parameter, 


le rs e080, 
Le 
constant = 3 


dimensionless ratio of 
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Fig. 1 Circular fin for radiant heat transfer 


of the fin.? The purpose of this paper is to solve the equations 
governing the temperature distribution and fin efficiencies for a 
nonlinear heat 


transfer, i.e., fourth-power law of absolute tem- 


peratures, and to present curves relating the fin efficiency to the 
dimensionless independent variables. In addition, a comparison 
of the fin efficiencies for convective linear heat-transfer data [6 
as shown by Kays and London iY} to those calculated for non- 


linear radiation heat transfer is given 


Treatment of the Problem 


The diffusion of heat in an annular fin of constant thickness re- 
duces to the problem of obtaining the steady-state temperature 
distribution in a finite hollow circular cylinder representing the 
annular fin. Several restrictions and assumptions are inherent 
in solving this fin problem 


1 The heat flow and temperature distribution throughout the 
fin are independent of time, i.e., the heat flow is steady state 
2 There is no heat source in the fin 


3 The 


t The fin thickness is so small compared with its radius that 


temperature at the base of the fin is uniform 


the temperature gradient across the width may be neglected, i.e., 


heat flow occurs in the radial direction only 


5 The fin radiates to a vacuum space which is at a constant 


arbitrary temperature of absolute zero 


6 The fin material is isotropic 


These restrictions combine with the phy sical laws of heat trans- 


fer to establish a second order nonlinear differential equation 


which describes the heat flow 


within the fin and a set of two 
boundary conditions, one each at the radial extremities of the fin 


The equations in dimensionless form are 


de l dé’ 


vyh'4 = 0 
dr ‘ 


* The peg is considered to be of infinite thermal conductivity; radia- 
tion from the exposed end of the peg is not considered in the solution 


presented 
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By definition the fin efficiency 7 equals Q/Qmax; therefore fin 
efficiency due to radiation, based on the solution for temperature 


from equations (1), is* 
l ; 
of, 0’*[r'(p 
n = 
p 


The solution to equations (1) is not known to be tabulated, 


— 1) + 1]dr’ 


+ 1 


and a numerical computation of it for different parametric values 
of p and y was necessary. The solution was accomplished with 
the use of an IBM 704 digital computer, and the adiabatic bound- 
ary condition atr’ = 5. equations (1), was checked zero to w ithin 
five significant figures. 

After 0’ as a function r’ is known, there are several methods of 
numerical integration that may be applied to find the fin efficiency 
One method of calculating the fin efficiency was derived explicitly 
by breaking the fin into a series of concentric annuli of annular 
distance Ar, and specifying a fin temperature @,’ at r,’ 


nAr : 

By definition the fin efficiency is thus given by 
Q 

~ Quas 


mweoO4 (r.2 — rr? 


n 


Qmax _ 


n=1 


In reduced form, equations (3 become 


I 
v+it (ov 
> r (4) 
Wp? 

In the limit as m approaches infinity, the dependence of the 
right-hand side of equation (4) on W disappears, and the right- 
hand side, now equalling 7, is a function of pand y¥ only, recalling 
that @’ is a function of p and y. Thus equation (4), shows that, 
within the degree of approximation associated with the selection 
of the value of m, 


n~Sfil¥, p, Y) 
n = fp, Y) 
fA, p, Y) = & flo, Y) 


Another quadrature method of approximating » was employed, 
that of Lowan, Davids, and Levenson [10]. This method uses 
the Gauss mechanical quadrature method with Legendre poly- 
nominal representation of the integrand. The two methods of 
numerical integration agreed to within 1 per cent, and: this tends 
to confirm the approximation in equations (5). While the Gaus- 
sian method is probably the more accurate, the summation 
method was used because of the simplicity of programming for 
the digital computer. 


Discussion 


Fig. 2 shows the fin efficiency as a function of two independent 


parameters, pand y. Thus if the two independent conditions of 


geometry and operating material constants are known, the fin 


efficiency can be determined from these curves. An interesting 


T 7 
‘Since r’ = 


*, the fin efficiency applies only to the fin 
jo = F 


and 


not to the exposed end of the peg. 
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Fig. 2 Comparison of convection and radiation fin efficiencies 


comparison of the work reported here with estimates of the fin 
efficiency based on linearizing the radiation rate, i.e., by replacing 
the fourth-power law by an equivalent convection constant, h’ = 
eo4;*, multiplied by the first power of the absolute temperature 


» 


and then using the available data cited [6, 9], is shown in Fig. 2. 


From the curves it can be seen that an estimate based upon 
linearizing the radiation rate will lead to a significant error unless 
V ¥ is less than 0.2. An attempt to adjust the estimate by using 
a mean or average temperature of the fin surface will lead to an 
even greater error since an average temperature will lead to a 
smaller y and thus to a higher efficiency instead of the correct 
one It therefore appears that there is no accurate way to de- 
termine the correct radiation rate except to « msider the fourth- 
power law of the absolute temperature 

Because of IBM 
work presented in this paper did not consider the case of a straight 
fin, p = 1.000 1.001 
and p = was found to be 
negligible 


recoding complexities in the program, the 
However, values were computed for p= 
1.1 and the difference in fin efficiency 

It is therefore reasonable to consider that a straight 


fin may be approximated by p 1.001 


Conclusions 


The one-dimensional steady-state solution for radiation from 
one side of an annular fin supplies needed design information for 


satellite thermal-power dissipating suriaces It is pre sented here 


in graphical form with the fin efficiency plotted versus the design 


eo6, 


mairameters (r r 
pa \ F i 


linearizing the radiation rate and using the 


and p The previous method ol 


convective analysis 


[6] introduces a significant error which can approach a maximum 
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of 60 per cent The results presented here permit correct compu- 
tation of the heat-transfer rates for radiation-controlled heat- 


dissipating fins of circular geometry. 
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be used in more elaborate formulas. The solidus should be used instead of the 
horizontal! line for fractions wherever possible. 


e A normal paper should not exceed 8 pages in the Journal. It is about 8000 
words or 20 double-spaced typewritten pages exclusive of drawings. A Tech- 
nical Brief should not exceed 750 words or about one page in the Journal inclusive 
of drawings. 


e Originals and two copies of figures must accompany the manuscript. Line 
drawings should not be larger than 8' . « 11 in. and should be planned for reduc- 
tion to column width. Lettering should be large enough to be clearly legible 
when the illustration is reduced. The originals of line drawings: must be in 
India ink on white or pale blue tracing paper or tracing cloth. Photographs of 
equipment or test specimens must be glossy prints and should be used sparingly. 
Captions for figures should be typed double-spaced and included as the last page 
of the manuscript. The figure number and author’s name should be written in 
the margin or on the back of each illustration. 

e Titles of papers should be brief. 

e Authors can obtain copies of the ASME Manual MS-4, “An ASME Paper,” 
from the Editor and are urged to do so before drafting their papers in final form. 
e Papers published in Journals of the Transactions of the ASME must be pre: 
sented at a Meeting, either in person or by title. 

e An author is entitled to 25 preprints free of charge (in the case of two authors, 
15 each; three or more authors, 10 each). Larger quantities of preprints or 
reprints can be ordered from Editorial Department, The American Society of 
Mechanical Engineers, 29 West 39th Street, New York 18, N. Y. Quotations 
will be sent on request. 


ASME Publications Committee 






































° * * - - 
+ » - «™ e - coal o eee 
n in “ 
- P FS ee 7 \, 
} ee] . F 1 ae 
op | -“. . F t 
- - es del Pr ‘ Sn . 
: s =T 
ve 5 ESE ies i re F s i Pa ‘ % 
. . 7 _ Pd ci 3 oad a 
nd pss oe oe m= . : 
oy ‘ . ; . ? ay 
z : & i 7 
oo ak . a} . F 
t~ ‘ i y 
™~e, a - : i 
~ % . . 
. * . 
-~ , 7 » , 
‘ 
5 Se j ~ . . 
; 7 6 ; 
- - ‘ 
* - 
‘ 
« . = —— i : : 
7 ; : cae . : 
Land 7 ¢ ' 
— ““ 
. ¢ 
* m a ; < 
Se. in. F 
Saas j ‘ * i Sua ® : * 
: § ah “~ i ¥ 4 ae = 
# t at ; iy : ae 
f x mcnamint fount 





